


P.BEAUFILS /M. LAMARCHE / Y.MUGGIANU 


PROGRAMMES DE 
MATHEMATIQUES SUR 


AMSTRAD 








EYROLLES 


PROGRAMMES DE 


MATHEMATIQUES SUR 
AMSTRAD 


CHEZ LE MEME EDITEUR 


Des mêmes auteurs : 
— Programmes de physique sur Amstrad. 


Autres ouvrages : 
C. DELANNOY  — Faites vos jeux avec Amstrad. 


— Je débute en BASIC Amstrad. 


P. BIHAN — Programmation sur AMSTRAD PCW 8256 et 8512 


Basic et fichiers. 
M. ROUSSELET — Ca/cul numérique sur Amstrad. 
©. LEPAPE L'assembleur facile du Z 80. 
Ph. DAX CP/M et sa famille. Guide d'utilisation. 


A. REVERCHON, M. DUCAMP 


— Mathématiques sur micro-ordinateur. 
1 - Analyse. 
2 - Algèbre. 


F. MANCHON, J.-M. NASR 


— Physique sur micro-ordinateur. 


J.-P. DELAHAYE 
— Dessins géométriques et artistiques avec votre micro- 


ordinateur. 


Nouveaux dessins géométriques et artistiques avec 
votre micro-ordinateur. 


M. JAMES Introduction à l'intelligence artificielle sur micro- 
ordinateur. 





PROGRAMMES DE 
MATHEMATIQUES SUR 


AMSTRAD 


par 


P.BEAUFILS / M. LAMARCHE / Y.MUGGIANU 





61, boulevard Saint-Germain — 75005 PARIS 
EYROLLES 1986 


Si vous désirez être tenu au courant de nos 
publications, il vous suffit d'adresser votre carte 
de visite au : 


Service « Presse », Editions EYROLLES, 
61, Boulevard Saint-Germain 
75240 PARIS CEDEX 05, 


en précisant les domaines qui vous intéressent. 
Vous recevrez régulièrement un avis de paru- 
tion des nouveautés en vente chez votre libraire 
habituel. 





« La loi du 11 mars 1957 n’autorisant, aux termes des alinéas 2 et 3 de l’article 41, d’une 
part, que les « copies ou reproductions strictement réservées à l’usage privé du copiste et non 
destinées à une utilisation collective » et, d’autre part, que les analyses et les courtes citations 
dans un but d’exemple et d'illustration, « toute représentation ou reproduction intégrale, ou par- 
tielle, faite sans le consentement de l’auteur ou de ses ayants droit ou ayants cause, est illicite » 
(alinéa 1° de l’article 40). » 


« Cette représentation ou reproduction, par quelque procédé que ce soit, constituerait une 
contrefaçon sanctionnée par les articles 425 et suivants du Code pénal. ». 


© Editions EYROLLES, 1986 


Préface 


Nous sommes physiciens et pourtant le titre de cet ouvrage parle de 
Mathématiques, ceci mérite une explication. 

Ce livre est initialement écrit par des physiciens pour des physiciens : 
— pour éviter qu'ils soient bloqués dans la résolution de leurs problèmes 
par des difficultés mathématiques et par là leur suggérer de s'attaquer à 
des problèmes encore plus complexes ; 
— pour leur fournir des méthodes de calcul numérique ou d'approxima- 
tions rapides et efficaces; 
— pour leur faire découvrir par l'exemple des domaines des mathémati- 
ques qu'ils connaissent peut-être mal. 


Certes les résultats fournis par l'ordinateur ne sont qu'approchés, 
mais le physicien préfère souvent une solution approximative plutôt 
qu'une absence de solution exacte. Il faut dire que les approximations 
faites dans les calculs sont parfois bien inférieures à celles consistant à 
négliger certains phénomènes physiques secondaires. 


Mais nous avouons avoir poursuivi également un autre but : intéres- 
ser nos collègues mathématiciens aux possibilités de l'ordinateur. Nous 
les soupçonnons en effet de penser que l'usage de l'ordinateur «tue» 
les mathématiques. Nous n'en croyons rien : la conception des program- 
mes et la critique nécessaire des résultats qu'ils fournissent ouvrent aux 
mathématiciens un grand champ de travail. 


Certes, nous concevons qu'ils bondissent lorsque l AMSTRAD appelle 
PI un nombre rationnel. Remarquez toutefois qu’il n'accepte pas PI ou 
4/2 à une demande d'introduction de données ou qu'il ne fournit pas PI 
pour résultat. Dans la mesure où l’on aboutira à une solution numérique 
approchée du problème, cet inconvénient n'est pas grave, mis à part la 
question de langage. 


Mais par d'autres côtés, l'utilisation de l'ordinateur fait progresser 
dans la précision de l'écriture. Ainsi sin? x + 3 doit disparaître au profit 
par exemple de (sin (x}) ? 2 + 3, ce qui est un bien. De même, l’usage 
du symbole * pour la multiplication évite des erreurs et permet l'emploi 
de variables désignées par deux lettres telles que «ab», ce qui est bien 
pratique. 


Quel peut être l'intérêt pour le mathématicien d’une solution numé- 
rique approchée ? Elle lui permet de déceler une solution sans intérêt 
ou par son invraisemblance de trahir une erreur de calcul. La solution 
numérique, le graphe d'une fonction par exemple, sont à l'expression 
littérale ce que la photographie d'identité est au visage : ce n’est pas le 
visage, mais c'est très utile pour reconnaître un visage. 


Bien qu'il soit spécialisé dans l'application numérique, l'ordinateur 
privilégie la notion de fonction, ce qui en fait en ce domaine un instru- 
ment pédagogique puissant. En voici un exemple : soit à dessiner sur 
l'écran des vaguelettes, ou un trèfle à quatre feuilles, ou un serpent ; si 
vous définissez ces figures point par point, il vous faudra introduire des 
centaines de nombres. Si vous les concevez comme graphes de fonction, 
il suffit d'écrire leur équation et de préciser le domaine de la variable. 
Toute personne ayant fait cette démarche sera convaincue de l'intérêt 
de la notion de fonction. 


Venons-en enfin à la critique des résultats. Si vous tracez le graphe 
de y = sin (4.x) en prenant pour x des valeurs entières ou si vous inté- 
grez cette fonction avec dx = 1, vous obtiendrez un résultat fantaisiste. 
Il est très difficile de réaliser un programme incollable et si le bon sens 
permet de rejeter la solution improbable, il faut parfois un sens criti- 
que et mathématique très aigu pour rejeter une solution apparemment 
convenable, mais en réalité fausse. 
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Ces programmes ont été écrits pour les ordinateurs 
AMSTRAD. Les commentaires qui les accompagnent 
permettent de les adapter facilement à d'autres ordina- 
teurs travaillant en BASIC. 


Les illustrations ont été réalisées avec l'imprimante 
AMSTRAD CMP 2000. 


La compréhension des commentaires requiert le 
niveau mathématique d'une classe terminale scientifique. 


Graphes de fonctions 
Y (X)ET R (A) 


Il ne s’agit pas ici de découvrir une fonction quelconque sur l'inter- 
Valle — © , + , mais d'obtenir une représentation graphique d'une 
fonction connue sur un intervalle bien déterminé. En particulier, la 
fonction doit être définie sur cet intervalle et ne pas présenter d'accident 
sur un domaine égal à 1/600 du domaine d'étude. Nous vous présentons 
donc un programme fait pour une fonction particulière, mais les com- 
mentaires qui suivent vous permettront de l'adapter facilement à toute 
autre fonction de votre choix. 


Tracé de Y, (X) 
Nous avons d'abord choisi de représenter la fonction numéro : 


ÿ1 =8.exp(—x).sin (8.x) (ligne 40) 
x variant de xa = 1 à xb = 3 (ligne 60) 


L'étendue /x du domaine d'étude fixe l'échelle en abscisse ex 
(ligne 120). L'intervalle entre deux valeurs de x considérées est dx 
(ligne 1030) = ;x/600. Le nombre 600 est fixé par le nombre de pixels 
sur une largeur d'écran compte tenu des marges. Il est inutile de calculer 
plus de valeurs car les points correspondants ne se distingueraient pas 
sur l'écran. 


aintenant le problème de l'échelle en ordonnées : 
Abordons maintenant | obl de l'échell d ées 


Le sinus est compris entre — 1 et + 1, l'exponentielle décroissante 
prendra sa valeur maximale pour x = 1, elle est inférieure à 0,37, multi- 
plié par 8, cela fait moins de 3. y est donc compris entre — 3 et + 3 
dans le domaine considéré. On peut bien sûr demander à l'ordinateur de 
déterminer ce domaine de valeur de y en calculant y pour les 600 valeurs 
de x et en recherchant les valeurs maximale et minimale. Ce procédé est 


utilisé dans les programmes «valeur moyenne» et «primitives», nous ne 
l'aborderons donc pas ici, mais remarquez que notre petit raisonnement 
mathématique a évité un calcul laborieux à l'ordinateur. 


Connaissant les valeurs extrémales yp et ym (ligne 80), on en déduit 
l'ampleur ;y de l'intervalle de valeur de y et compte tenu de la hau- 
teur de l'écran (360 en ligne 150), l'échelle ey . [[ by est le décalage qui 
va nous permettre de placer l'axe des abscisses. Si la fonction étudiée 
garde le même signe, positif, par exemple, on pourra choisir ym = 0 
afin d'observer l'axe des abscisses sur l'écran. 

Tout ceci étant fait, la fonction est calculée pour chaque valeur de 
x (ligne 1050). Les coordonnées du point représentatif sont calculées 
compte tenu des marches, de l'échelle et du décalage en ligne 1080 et 
1090 (fig. 1). 





1 pixel 
———————————— 


600 pixels 








Fig. 1. 


Si vous supprimez les lignes 1302 à 1306, vous obtenez sur l'écran 
une distribution de points suggérant la courbe. En ajoutant ces lignes, 
vous joignez ces points par des droites (théoriquement, en fait ce sont 
des courbes en escalier). À vous de choisir entre une suggestion de cour- 
be assez exacte et une courbe imposée un peu moins exacte (fig. 2). 


Fig. 2. 


Tracé de Y, (X) 


Lorsque la flèche apparaît, vous pouvez appuyer sur la touche cor- 
respondante pour observer le graphe de la deuxième fonction obtenu 
de la même manière que le précédent (fig. 3) : 

y = exp (x/8) . sin (2 .x) 


x variant de —3 à +3 





Fig. 3. 


Tracé de Y: (Y, (X)) 


L'étape suivante du programme montre (fig. 4) le graphe de la fonc- 
tion composée à partir des deux précédentes, x variant de 1 à 3. Lorsque 
la variable prend la valeur x, la première fonction en donne y, qui est la 
valeur de la variable pour la seconde fonction (ligne 260). 


Fig. 4. 


Tracé d’une courbe définie paramétriquement 


Un point de cette courbe a pour abscisse y, (x) et pour ordonnée 
y: (x), nous avons représenté figure 5 la portion de courbe correspon- 
dant à une variation de x de 1 à 3. 


Structure du programme 


Lignes) 


40 
60 
170 
1000 
200 
340 


définition des fonctions 

préparation des graphes 

envoi à 1000 pour le tracé de la courbe 
tracé de y, ou de y; 

tracé de la fonction composée 

tracé de la courbe paramétrique 


Fig. 5. 


Liste des variables 


x variable 

n numéro de la fonction 

y valeur de la fonction 

xa première valeur de x 

xb dernière valeur de x 

ym, yp Valeurs minimale et maximale de y 

ex,ey échelles en x et en y 

ix, iy  largeurs des domaines de valeurs de x et y 
bx,by décalages en abscisses et en ordonnées 

dx _ intervalle entre deux valeurs de x consécutives 
xe,ye coordonnées du point sur l'écran 

XxXx,yy coordonnées du point précédemment marqué 
v valeur de la variable dans le cas des fonctions composées 
dv écart entre deux valeurs de v 


Premier Programme 


10 REM'"graphes de fonctions" 
20 REM 

SO CLS:GOTO 60 

40 y(1)=8*EXP (-x) #SIN(B#x) 
45 RETURN 

90 y(2)=EXP (x/8) #xSIN(2#x) 

55 RETURN 

60 xa(1)=1:xb(1)=3 

70 xa(2)=-3:xb(2)=3 

80 ym(1)=-3:yp(1)=3 

90 ym(2)=-1.5:yp(2)=1.5 

100 FOR n=1 TO 2 

110 ix(n)=xb(n)-xa(n) 

120 ex (n)=600/ix(n) 

130 bx(n)=-xa (n)#ex (n) 

140 iy(n)=yp{(n)-ym(n) 

150 ey(n)=360/iy(n) 

160 by(n)=-ym(n)#ey (n) 

170 GOSUB 1000 

180 NEXT n 

190 REM 

200 MOVE 0,20+by(2):DRAWR 640,0 
210 MOVE 20,0:DRAWR 0,400 
220 dv=ix(1)/600 

230 FOR v=xa(1) TO xb(1) STEP dv 
240 x=vi:xe=20+bx (1) +ex (1) #x 
250 GOSUB 40 

260 x=y(1) 

270 GOSUB 50 

280 ye=20+by (2) +ey (2) #y (2) 
290 PLOT xe,ye 

292 IF v=xa(i) THEN GOTO 296 
294 DRAWR xx-xe,yy-ye 


296 
300 
310 
329 
330 
340 
350 
360 
580 
390 
409 
410 
420 
430 
432 
434 
436 
440 
970 
980 
990 
1000 
1019 
1020 
1030 
1040 
1050 
1960 
1070 
1080 
1090 
1300 
1302 
1304 
1306 
1310 
1320 
1330 
1340 
1350 


xx=xe: yy=ye 
NEXT v 
LOCATE 39,25:PRINT CHR#$(243) 
q=INKEY(1):1F q<>0 THEN 320 
REM 
CLS 
MOVE 320,0:DRAWR 0,400 
MOVE 0,200:DRAWR 640,0 
FOR x=xa(1) TO xb(1) STEP ix(1)/600 
GOSUB 40 
GOSUB 50 
xe=320+tey(1)#y(1) 
ye=200+ey (2) #y(2) 
PLOT xe,ye 
1F x=xa(1i) THEN GOTO 436 
DRAWR xx-xe,yy-ye 
xx=xe:yy=ye 
NEXT x 
END 
REM 
REM 
MOVE O0,20+by(n):DRAWR 640,0 
MOVE 20,0:DRAWR 0,400 
LOCATE 35,24-(by(n)/16)+1:PRINT STR#&(xa(n)) 
LOCATE 36,24-(by(N)/16)+1: PRINT STR#&(xb (n)) 
LOCATE 1,24-(by(n)/16)-1:PRINT"O" 
dx (n)=ix (n) /600 
FOR x=xa(n) TO xbin) STEP dx(n) 
GOSUB 40:GOSUB 50 
ke=20+bx (n) +ex (n) #x 
ye=20+by(n)+ey(n)#y(n) 
PLOT xe,ye 
IF x=xa(n) THEN GOTO 1306 
DRAUR xx-xe,yy-ye 
xx=xe:yy=ye 
NEXT x 
LOCATE 39,25:PRINT CHR$ (245) 
q=INKEY(1):1F q<>0 THEN 1330 
CLS 
RETURN 


Coordonnées polaires 


La courbe est tracée après avoir choisi l'échelle (e = 50 en ligne 80) 
et le domaine de variation de l'angle (ligne 60 à 70). Le point est placé 
en lignes 130 à 150. Il est relié au précédent en lignes 160 à 180 (fig. 6). 


Vous prendrez garde que certaines imprimantes, dont celle-ci, défor- 
t l’image. Pour obtenir sur le papier une image non déformée, il 
multiplier ye par un coefficient adéquat. Vous mesurerez le coeffi- 
cient de distorsion en faisant recopier par l'imprimante un carré tracé 


men 
faut 


sur | 


‘écran. 


Fig. 6. 


Deuxième Programme 


10 REM"COURBE EN POLAIRE" 
20 REM 

30 CLS:RAD:GOTO 60 

40 R=4#SIN(a/2)+SIN(2%a) 
50 RETURN 

60 ai=-2*PI 

70 a2=2%PI 

80 e=50 

90 ia=a2-ali:da=ia/500 

100 MOVE 320,200:DRAWR 320,0 
110 FOR a=ai TO a2 STEP da 
120 GOSUB 40 

130 xe=320+e#r#C0S (a) 

140 ye=200+e#r#SIN(a) 

150 PLOT xe,ye 

160 1F a=al THEN GOTO 180 
170 DRAWR xx-xe,yy-ye 

180 xx=xe: yy=ye 

190 NEXT a 


Dans la courbe de la figure 6, on peut reconnaître un demi cœur 
correspondant à une variation de a de — pi à pi. Tournons la courbe de 
pi/2, modifions le domaine d'étude en conséquence, prenons la valeur 
absolue de l'angle pour obtenir la courbe symétrique et nous avons le 
cœur. En fait nous avons été amenés à modifier quelque peu l'équation 
pour creuser le cœur dans sa partie inférieure (fig. 7). Les corrections 
à apporter sont donc les suivantes : 


150 


190 


Fig. 7. 


Troisième Programme 


REM"COURBE EN POLAIRE" 

REM 

CLS:RAD:GOTO 40 

aa=ABS (a-P1/2) 

R=4XSIN(aa/2) +SIN(2#aa)#(1+0.3#SIN(aa) ) 
RETURN 

ai=-P1/2 

a2=P1+#3/2 

e=50 

ia=a2-al:da=ia/500 


MOVE 320,200:DRAWR 320,0 
FOR azai TO a2 STEP da 
GOSUB 40 
xe=320+e*#r*C0S (a) 
ye=200+e*r#+SIN(a) 


PLOT xe,ye 

IF a=ai THEN GOTO 180 
DRAWR xx-xe,yy-ye 
xx=xe: yy=ye 

NEXT a 


Graphisme mathématique 


LES DIFFÉRENTS SYSTÈMES DE COORDONNÉES 


Il y a trois principaux systèmes de coordonnées que l'on emploie 
pour représenter graphiquement les fonctions et les courbes. 


1. Les coordonnées cartésiennes : 


Un point M du plan est déterminé par son abscisse X et par son 
ordonnée Y. 
On posera OH = X et HM = Y. 


Y 





L'équation est alors donnée soit sous forme explicite Y = f (X), 
soit sous forme implicite f (X, Y) = O. 


2. Les coordonnées polaires : 


Les coordonnées polaires d'un point M sont définies à partir d'un 
axe OX, appelé axe polaire origine, contenant le point O appelé pôle. 


Si OY est un axe passant par M, on définit l'angle des axes : 
(OX, OY) = « et la mesure algébrique OM = p. 


a et p constituent les coordonnées polaires du point M. 
La relation p = f (a) est l'équation polaire d’une courbe. 


3. Les coordonnées paramétriques : 


La position d'un point M est donnée par son abscisse X et son 
ordonnée Ÿ, mais chacune d'elles étant exprimée en fonction d'un para- 


mètre, lequel est en général un angle «, ou le temps t. 


Programme 1 


Ce programme va vous permettre de dessiner la courbe de toutes les 


fonctions en coordonnées cartésiennes. 


Déroulement du programme 


Lignes) 
40 l'ordinateur demande une chaîne AS$ qui est la fonction 
à étudier 
50 l'ordinateur vous demande de lui fournir un nombre N. 


Il tracera la courbe pour les valeurs de la variable compri- 


ses entre — N et + N (voir ligne 100) 


80 à 200 l'ordinateur ouvre un tableau de 640 valeurs. La fonction 
est donnée en ligne 120. Une échelle automatique est 
calculée par l'ordinateur (ligne 130). En ligne 110, on 
rentre les conditions particulières nécessaires à l'étude 


de la fonciton 
300 à 330 on trace les axes 


340 à 430 l'ordinateur trace la courbe en continu (et non point par 


point) 


REM FONCTION Y=F(X) 
REM 66 DEEE DE DE DE DE JE IE JE 
MODE 2 

INPUT "ETUDE DE LA FONCTION Y=";A$ 
INPUT "INTERVALLE D'ETUDE N=";N 
CLS 
LOCATE 1i5,12:PRINT "CALCUL EN COURS" 
DIM A(640):S=0 
FOR F=0 TO 640 
LOCATE S,1S5:PRINT F 

X=(F-320)#N/320 

IF X=0 THEN GOTO 200 
A(F)=X#SIN(1/X) 

IF ABS(A(F))5S THEN S=AES (A(F)) 
NEXT F 

CLS 

REM TRACE DE LA FONCTION 

REM 36606266 DE DE DE DE DE DE DE JE DE DE DE DH 

FLOT 520,0: DRAWR 0,400 

PLOT 0,200:DRAWR 640,0 

T=0 

FOR F=0 TO 640 

IF AES(A(F))/S=1 THEN T=d:GOTO 400 
Y=A(F) /S#190 

IF T=0 THEN FLOT F,Y+200:T=1:GOTO 400 
DRAWR 1,Y-—AY 

AY=CINT(A(F)/S#190) 

NEXT F 

LOCATE S,22:FRINT "Y=";:A#$ 


LOCATE 3,23:PRINT "INTERVALLE D'ETUDE :"3;3-N;" 





Y=X2SIN(1/X) 

INTERVALLE D’ETUPE :-1 , 1 
Ready 
CALL &48688 


11 


12 


Y=SIN(X) 
pe DTERVALLE D'ETUE 1-6.28,, 6.20 
CALL 40088 7 





Y=LOC(ABS(A2+5)) 
INTERVALLE D'ÉTUDE :-5 , S 
Ready 
CALL 44888 


des 


Programme 2 


Il présente l'étude d'une fonction Y = TAN (X) qui comprend donc 
branches asymptotiques. 


Observez particulièrement les lignes 125 et 360. 


RKEM FONCTION Y=F(X) 
REM #36 3636 EEE EE DE DE DEEE 


MODE 2 
INFUT "ETUDE DE LA FONCTION Y=";A% 
INFUT "INTERVALLE D'ETUDE N=";N 
CLS 
LOCATE 15,12:PRINT "CALCUL EN COURS" 
DIM A(640) : S=0 
FOR F=0 TO 640 
LOCATE 3,15:PRINT F 
X= (F-35320) *N/320 
A(F)=TAN(X) 
IF ABS(A(F))>50 THEN A(F)=0:GOTO 200 


IF ABS(A(F))>S THEN S=ABS(A(F)) 
NEXT F 

CLS 

REM TRACE DE LA FONCTION 


REM +36 36262666 DE DE DE DE ED DE EE 

FLOT 320,0:DRAWR ©,400 

FLOT 0,200:DRAWR 640,0 

T=0 

FOR F=0 TO 640 

IF A(F)=0 OR ABS(A(F))/S=1 THEN T=0:GOTO 400 
Y=A(F) /S*190 

IF T=0 THEN FLOT F,Y+200:T=1:GOTO 400 

DRAWR 1,Y-AY 

AY=CINT (A(F)/S*#190) 

NEXT F 

LOCATE 3,22:FRINT "Y=":A# 

LOCATE 3,23:FRINT "INTERVALLE D'ETUDE :";-N;", 





D'EIUDE :-6.28 , 6.28 
Ready Q 
CALL 44088 
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Y=TAMCXD/X 
INTERVALLE D'’ETUDE :-S, 
HT 





Programme 3 


Ce programme vous permettra de calculer l'intersection de la courbe 
avec l'axe des X par la méthode de dichotomie. 


MÉTHODE DE DICHOTOMIE 


Considérons l'équation sin X = 0.2 X — 0.5, la variable X étant 
exprimée en radians. On ne sait pas la résoudre par les méthodes usuelles 
algébriques. Habituellement on trace Y = sin X, Z = 0.2 X —0.5et on 
a la solution graphique au point d'intersection. 





Fig. 1. 


Une autre méthode consiste à essayer des valeurs et à tenter de 
déterminer pour quelle valeur de X, sinX s'approche de 0.2 X — 0.5, ou 
encore pour quelle valeur de X l'expression : 

«sin X — 0.2 X + 0.5» est nulle. 


Pour certaines valeurs de X, l'expression précédente (appelons-la 
f (X)) est positive, pour d'autres elle est négative. On a ainsi deux 
bornes de X, soit À pour l’inférieure et B pour la supérieure entre les- 
quelles existe la racine. L'astuce consiste à resserrer la fourchette (A, B). 
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Remarque importante : 


La racine se trouve dans l'intervalle À, B si f (A) et f (B) sont de 
signes contraires ; ou encore si le produit f (A).f (B) est négatif. 


La dichotomie permet de réduire la fourchette jusqu'à ce que 
|B — A] soit aussi petit que l’on veut. On écrit |B — A| < PR(précision). 









A 


borne 
inférieure 


B 


borne 
supérieure 


Fig. 2. 


A+B 
Procédé : Calculons 2 = 








Fig. 3. 


A+B 
— Si f (==) .f (A) < 0, cas 
de la figure ci-dessus, la solution se 


B 
. Les bor- 





trouve entre A et 


nes sont maintenant À comme in- 





férieure et comme supé- 
rieure. Au fond, si f (X).f (A) <0 
avec X = 


B par X. 





, il fautremplacer 


nouvelle borne. Deux cas sont possibles. 








Fig. 4. 
| A +B 
ss ).ftA)>0 te 
solution se trouve entre la nouvelle 
+ B 
borne inférieure et l’an- 





cienne borne supérieure B. Si tra- 
ditionnellement la borne inférieure 
est désignée par «A», on peut dire 


+ B 
que le nouveau «A » est 


Au fond, si f (X).f (A) > 0 avec 


A+B . 
= Dur il faut remplacer A 





par X. 


— On recommencera l'opération jusqu’à obtenir la précision voulue. 


Déroulement du programme 


Quelques variantes par rapport aux programmes précédents : 


Ligne(s) 


55 l'ordinateur pose la question. concernant la graduation 
des axes (OUI ou NON) 
520 à 560 si B$ — "OUI ”, on gradue les axes au 1/5 et au 3/5 de 
N (ligne 56) (Procédure un peu longue) 
1000 à 1110 méthode de dichotomie. PR : précision prise ici (au 
1/1000) 


Remarque importante 
Dans certains cas, il peut se faire que Q = FNY (X) (1060) soit 
égal à zéro. Rajouter alors en 1065 : 


Ligne 
1065 IF Q = O0 THEN PRINT "X ="; X : STOP 


10 REM FONCTION Y=F(X) 
20 REM 434646062606 36226 De 


30 MODE 2 

40 INPUT “ETUDE DE LA FONCTION Y=";A$ 
30 INPUT "INTERVALLE D'ETUDE N='";N 

95 INPUT "GRADUATION DES AXES";B#$ 

56 IF B$="OUI" THEN P=N/5 

60 CLS 

70 LOCATE 15,12:PRINT "CALCUL EN COURS" 
80 DIM A(640):S=0 

90 FOR F=1 TO 640 

95 LOCATE 3,15:PRINT F 

100 X=(F-320) #N/320 

120 DEF FNY(X)=SIN(X)-0.2#X+0.5 

125 A(F)=FNY(X) 

1350 IF ABS(A(F))2S THEN S=ABS(A(F)) 

200 NEXT F 

210 CLS 

300 REM TRACE DE LA FONCTION 

S10 REM 446262626266 0606 ED DE UUE DE DD DE DE 

320 PLOT 320,0: DRAWR 0,400 

330 PLOT 0,200:DRAWR 640,0 

340 T=O:1F B$="OUI" THEN GOSUB 520 

350 FOR F=i TO 640 

S60 1F ABS(A(F))/S=1 THEN T=0: GOTO 400 
370 Y=A(F)/S#190 

380 1F T=O0O THEN PLOT F,Y+200:T=1:GOTO 400 
390 DRAWR 1,Y-AY 

400 AY=CINT(A(F)/5#190) 

410 NEXT F 

420 LOCATE S,20:PRINT "Y=";A$ 

430 LOCATE S,21:PRINT "INTERVALLE D'ETUDE :";-N;", 


450 LOCATE 3,3:PRINT “INTERSECTION AVEC L'AXE X'":1 
NPUT C$ 
460 IF C$="OUI" THEN GOTO 1020 ELSE STOP 


500 
510 
520 
530 
9540 
990 
560 
970 
1000 
1010 
1020 
1030 
1040 
1050 
1060 
1070 
1080 
1090 
1100 
1110 


REM GRADUATION DES AXES 
REM XD DEEE E DE DEEE NE DE 
FOR F=i TO 640 
FOR L=-N TO N STEP P 

IF F=CINT(640%xL/N) THEN PLOT F,200:DRAWR O,-5 
NEXT L 
NEXT F 
RETURN 

REM DICHOTOMIE 

REM #66 

INPUT "PRECISION PR=";PR 

INPUT "VALEUR INF. DE L'INTERVALLE A=";A 
INPUT "VALEUR SUP. DE L'INTERVALLE B=";B 
X=(A+B) /2 

O=FNY (X) : R=FNY (A) 

Z=G*R 

IF Z>20 THEN A=X 

1F Z<O THEN B=X 

IF ABS(B-A)<PR THEN PRINT “"X=";X:STOP 
GOTO 1050 










S HersECHION auRe L'AE x 


PRÉCISION PR:? .681 
VALEUR pe DE L’INTER 
VALEUR SUP. DE L’INTER 
X= .99902344 






Nrrensterion AVEC L'AXE X 
9 QUI 


RER TON PR=? .691 
VALEUR . D L'INTERVALLE A=? 2 
VALEUR UP: DE L'INTERVALLE B=? 3 


X= 2.99982344 
Break in 1184 

ady 
GOTO 1839 

LEUR INF. DE L’INTERVALLE A=? 8 
VALEUR SUP. DE L’ INTERVALLE B=? - 
X=-8.692382813 LR à 


YESIN(X) -8.24X+8.5 
INTERVALLE D'ETUDE :-3.5 , 3,5 
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Programme 4 


Ce programme va vous permettre de faire apparaître sur l'écran, 
si vous le désirez, la dérivée Y’ = f (X) de la fonction Y = f (X). 


Remarque importante : Dans les programmes 3 et 4, n'oubliez pas 
d'écrire en ligne 110 les conditions particulières à l'étude de la fonction 
(intervalle de définition, etc.). 


La dérivée est calculée numériquement à droite (ligne 140) et à 
gauche (ligne 150) de chaque point considéré et ensuite moyennée 
(ligne 160). 


Ligne (s) 
170 et 180 calcul de l'échelle automatique pour visualiser la courbe 
et sa dérivée 
440 l'ordinateur demande à l'utilisateur s’il veut visualiser 
la dérivée Y’ = f (X) 
450 l'ordinateur demande à l'utilisateur le calcul de l'inter- 


section avec l'axe des X (dichotomie) 


Les programmes 5, 6, 7, 8 et 9 présentent des études de courbes en 
coordonnées polaires. Les méthodes de calcul sont identiques aux 
programmes 1 et 2. 


10 REM FONCTION Y=F(X) 

20 REM #22 26062606 D ED DE 

SO MODE 2 

40 INPUT "ETUDE DE LA FONCTION Y=";A# 
90 INPUT "INTERVALLE D'ETUDE N='";N 

393 INPUT "GRADUATION DES AXES'";E$ 

56 IF B$="OUI" THEN P=N/5 

69 CLS 

70 LOCATE 15,12:PRINT "CALCUL EN COURS" 
80 DIM A(640):DIM D(640) :S=0:J=0: DX=0. 001 
90 FOK F=1 TO 640 

95 LOCATE 3,15:PRINT F 

100 X=(F-320)#N/320: A=X 

120 DEF FNY(X)=SIN(X)-0.2#X+0.5 

130 A(F)=FNY(X) 

140 X=A+D%:Z1=FNY(X) 

130 X=A-DX:Z2=FNY(X) 

160 D(F)=(Z1-22)/(2#DX) 

170 IF AHS(A(F))2>S THEN S=<ABS(A(F)) 
180 IF ABS(D(F))2>J THEN J=ABS(D(F)) 
200 NEXT F 

210 CLS 

300 RKEM TRACE DE LA FONCTION 

310 REM 6666 DEEE HE DE DE DE DE JE 

320 PLOT 320,0:DRAWR 0,400 

330 PLOT 0,200: DRAWR 640,0 

340 T=0:1IF B+="OUI" THEN GOSUR 520 
330 FOR F=i TO 640 
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360 
570 
580 
390 
400 
410 
20 
439 
440 
450 
455 
460 
5300 
510 
920 
930 
940 
990 
560 
370 
600 
610 
620 
650 
640 
650 
660 
670 
700 
710 
800 
1000 
1010 
1020 
100 
1040 
1050 
1060 
1070 
1080 
1090 
1100 
1110 


IF ABS(A(F))/S=1 THEN T=0:GOTO 400 

Y=A(F) /S5#190 

IF T=0 THEN PLOT F,Y+2R00:T=1:GOTO 400 

DRAWR 1,Y-AY 

AY=CINT(A(F)/S#190) 

NEXT F 

COCATE 3,20: PRINT "Y=":A#$ 

LOCATE 5,21:PRINT "INTERVALLE D'ETUDE :"3;3-Ns;","3N 
LOCATE 3,1:FRINT "VISUALISATION Y‘=F(X)":INFUT D# 
LOCATE 5,4:FRINT "INTERSECTION AVEC L'AXE X'":INFUT L 
IF D#="OUI" THEN GOTO 420 

IF C#="OUI" THEN GOTO 1020 ELSE STOF 
REM GRADUATION DES AXES 

REM HD DE DEMO DE HE ME DE DE DE DEEE 
FOK F=i TO 640 

FOR L=-N TO N STEF PF 

IF F=CINT(640#%L/N) THEN PLOT F,200:DRAWR O,-5 
NEXT L 

NEXT F 

RETURN 
KEM DERIVEE Y'=F(X) 

REM HUE 

T=0Q 

FOR F=1 TO 440 

IF AHS(D(F))/J=1 THEN T=0:GOTO 700 

D=D (F) /J#190 

IF T=0 THEN FLOT F,D+200:T=1:GOT0 700 

DRAWR 1,D-AD 
AD=CINT(D(F)/J#190) 

NEXT F 

IF C+$="NON" THEN STOP 

REM DICHOTOMIE 

REM #66 ee EE 


INPUT "PRECISION FR="3; FR 
INFUT "VALEUR INF. DE L'INTERVALLE A=";A 
INFUT "VALEUR SUF. DE L'INTERVALLE EB=";E 
X= (A+E) /2 

Q=FNY (X) : R=FNY (A) 

Z=CO*R 


IF Z:0 THEN A=X 

IF Z£O THEN EB=X 

IF ABS(B-A)<FR THEN FRINT "X=";X:STOF 
GOTO 1050 


VISUALISATION Y’=F(X) 
9 QUI 







a NTRRSECTION AVEC L’AXE X 
ISION PR=? 


PRÈC .B01 
VALEUR BUP. BE L'INTÉRUALÉE 82 
X=-1,99902344 






1 
è 











X)-Q<2#X+8.5 
D'ETUDE :-3 , 3 
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10 
29 
30 
40 
50 
99 
56 
69 
70 
80 
90 
100 
110 
120 
170 
140 
1350 
155 
160 
170 
175 
180 
190 
200 
210 
220 
270 


240 


REM COURBE : COORDONNEES FOLAIRES 
REM #26 2626363636 D DE DD DEEE DE DE DE ED DE DE DE DE DE DE HE 


REM 
CLS: MODE 2 

INPUT "ETUDE DE LA FONCTION: Z=";7# 
INFUT "VALEUR DE M=";M 

INFUT "VALEUR DE N=";N 
LOCATE 20,12:FRINT "CALCUL EN COURS" 
DIM A(BO00):DIM X(800):DIM Y (800) 
CX=0:CY=0:S=0 
FOR F=1 TO 800 

T=50%F /800 

LOCATE Z,1S:PRINT "T="3;CINT(T) 

DEF FNY(T)=SIN(M#T/N) Û 
ACFI=FNY CT) 2: X (F) =A(F) COS (TD: Y (F)=A CF) #SINCT) 
IF ABS(X(F))>CX THEN CX=ABS(X(F)) 

IF AERS(Y(F))>CY THEN CY=ABS (Y(F)) 
NEXT F 

IF CX>CY THEN S=CX 

IF CY>CX THEN S=CY 

CLS 

REM TRACE DE LA COURBE 

REM #66 6e DEEE EEE DE DE DE DE DE 

FLOT 320,0:DRAWR 0,400 

FLOT 0,200:DRAWR 640,0 

FOR F=1 TO 800 

IF F=800 THEN GOTO 400 

PLOT X(F)/S*180+320,Y(F) /S*#180+200:DRAWR (X(F+1)-X(F))/S*180, 


(Y(F+1)-Y(F))/S#180 


250 
400 
410 
429 
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NEXT F 

LOCATE 3,22:PRINT "Z=":2# 
LOCATE 3,23:FRINT "M=":M 
LOCATE 3,24:FRINT "N="3;N 










2 
TN 





10 
20 
30 
40 
30 
60 
70 
80 
90 
100 
110 
120 
150 
140 
150 
155 
160 
179 
175 
180 
190 
200 
210 
220 
250 
240 
1)- 
250 


400 
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REM COUREE : COORDONNEES POLAIRES 
REM 3636362636 46 66 DD DEEE DE DE DE DE DEEE DE DE DE DEEE NE 


REM 
CLS: MODE 2 

INFUT "ETUDE DE LA FONCTION:Z=";2# 
LOCATE 20,12:PRINT "CALCUL EN COURS" 
DIM A(90):DIM X(90): DIM Y (90) 
CX=0: CY=0: S=0 
FOR Fri TO 90 

T=5#F/90 

LOCATE 3,15:PRINT "T="3;CINT(T) 

DEF FNY(T)=SIN(T)#COS(2*#T) 

A CF) =FNY CT) 2 X CF) =A CF) #COS (TD 2 Y (FDA (FI XSINCT) 
IF ABS(X(F))>CX THEN CX=ABS(X(F)) 
IF AES(Y(F))>CY THEN CY=ABS(Y(F)) 
NEXT F 

IF CX>CY THEN S=CX 

IF CY)CX THEN S=CY 

CLS 

REM TRACE DE LA COURBE 

REM 46260606 D 6 6 Dé DE EME DE DE DE 

FLOT 220,0:DRAWR 0,400 

PLOT 0,200:DRAWR 640,0 

FOR F=1 TO 90 

IF F=90 THEN GOTO 400 

PLOT X(F)/S#180+320,Y(F)/S*180+200:DRAWR (CX(F+ 
X(F))/S#180, (Y(F+1)-Y(F))/S#180 
NEXT F 

LOCATE 3,22:PRINT "Z="3;7# 





Z=SIN(T)#COS(2#T) 


100 
110 
120 
130 
140 
150 
155 
160 
170 
175 
180 
190 
200 
219 
220 
230 
240 
1)- 
250 
400 
410 
429 


REM COURBE : COORDONNEES POLAIRES 
REM 626462626266 DE DE DE DEEE DEEE DE DE DE DD DE DE DE DE DEN 
REM 
CLS: MODE 2 

INPUT “ETUDE DE LA FONCTION:Z=";2# 
INPUT "VALEUR DE A=";A 

INPUT "VALEUR DE R="3;R 
LOCATE 20,12:PRINT "CALCUL EN COURS" 
DIM A(800):DIM X(800):DIM Y (800) 
CX=0:CY=0:S=0 
FOR F=1 TO 800 

T=2*P1+#F/800 

LOCATE 3,15:PRINT "T="3;CINT(T) 
DEF FNY(T)=A+2#R#C0S (T) 
A(F)=FNY(T):X(F)=A(F) COST): Y(F)=A(F)#SINUT) 
IF ABS(X(F))2CX THEN CX=ABS{(X(F)) 
IF ABS(Y(F))>CY THEN CY=ABS(Y(F)) 
NEXT F 

1F CX>CY THEN S=CX 

IF CY>CX THEN S=CY 

CLS 

REM TRACE DE LA COURBE 

REM 4636262666 DE DEEE DE DE DE DE DEEE 

FLOT 320,0: DRAWR 0,400 

PLOT 0,200:DRAWR 640,0 

FOR F=i TO 800 

IF F=800 THEN GOTO 400 

FLOT X(F)/S#180+320,Y(F)/S#180+200:DRAWR (X(F+ 
X(F))/S#180,(Y(F+1)-Y(F))/S#180 
NEXT F 

LOCATE 3,22:PRINT "Z=";2# 

LOCATE 3,23: PRINT "A=";A 

LOCATE 3,24: PRINT "R='";R 





Z:A42R#COS(T) 
R= i 
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Z=A+2#RACOS(T) 
a= 8 
R= 4 





Z=A+2#RACOS(T) 
A= 8 
RS 





Z=A+2#RaCOS(T) 
8= 18 
R= 1 


10 
20 
30 
40 
90 
93 
680 
70 
80 
90 
100 
110 
120 
130 
140 
150 
155 
160 
170 
175 
180 
190 
200 
210 
220 
230 
240 
1)- 
250 
400 
410 


REM COURBE : COORDONNEES FOLAIRES 
REM 3266 EEREOHEOHEHEONHE DH HE HE HEDE DE DE DE DE DE DEEE 
REM 
CLS: MODE 2 
INPUT "ETUDE DE LA FONCTION: Z=";2# 
INPUT "VALEUR DE R="3R 
LOCATE 20,12:FRINT "CALCUL EN COURS" 
DIM A(800):DIM X(800):DIM Y (800) 
CX=0:CY=0:S=0 
FOR F=1i TO 800 

=2*PI#F/800 

LOCATE 3,15:PRINT "T="3;CINT(T) 
DEF FNY(T)=R#T 
A(FI)=FNY(T):X(F)=A(F) #COS(T):Y(F)=ACF)#SIN(T) 
IF ABS(X(F))>;CX THEN CX=AES(X(F)) 
IF ABS(Y(F))2>:CY THEN CY=ABS(Y(F)) 
NEXT F 

IF CX>:CY THEN S=CX 

IF CY>CX THEN S=CY 

CLS 

REM TRACE DE LA COURBE 

REM HD 6 EME DE DE DE DE DE DE DE DE DE JE 

FLOT 320,0:DRAWR 0,400 

PLOT 0,200: DRAWR 640,0 

FOR F=1i TO 800 

IF F=800 THEN GOTO 400 

PLOT X(F)/S#180+320,Y(F)/S#180+200:DRAWR (X(F+ 
X(F))/S#180,(Y(F+1)-Y(F))/S#180 
NEXT F 

LOCATE 3,22:PRINT "Z="37# 

LOCATE 3,23:FRINT "R='"3R 
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10 
20 
30 
40 
90 
99 
56 
60 
70 
80 
90 


100 
105 
110 
120 
125 


176 


REM COUREE : COORDONNEES FOLAIRES 
REM 236626 DEEE DEEE DD DE DEEE DH DE DEEE DEN HN 
REM 
CLS: MODE 2 
INFUT "ETUDE DE LA FONCTION: Z="3:27# 
INFUT "VALEUR DE A=";:A 
INFUT "VALEUR DE EB=":HB 
LOCATE 20,12:FRINT "CALCUL EN COURS" 
DIM A(800):DIM X(800):DIM Y (800) 
CX=0:CY=0:S=0 
FOR F=i TO 800 
=2*P1I1*F/800 
IF COS(T)=0 THEN A(F)=0:GOTO 155 
LOCATE 3,1S:FRINT "T="; CINT(T) 
DEF FNY(T)=A+B/COS(T) 
ACF)=FNY(T) 
IF ABS(A(F)):3 THEN A(F)=0:X(F)=0:Y(F) =0:GO0OTO 





155 


130 


X(F)=A(F) #COS(T):Y(F)=A(F)#SIN(T) 


140 IF ABS(X(F)):CX THEN CX=ABS(X(F)) 
150 IF ABS(Y(F))2>CY THEN CY=ABS(Y(F)) 
155 NEXT F 


160 


IF CX:CY THEN S=CX 


170 IF CY>:CX THEN S=CY 

175 CLS 

180 REM TRACE DE LA COURBE 
190 REM #26 36 626 DEEE DEEE DEEE 

200 PLOT 320,0: DRAUR 0,400 

210 PLOT 0,200:DRAWR 640,0 

220 FOR F=1 TO 800 

225 IF A(F)=0 THEN GOTO 250 

230 1F F=800 THEN GOTO 400 

249 FLOT X(F)/S#180+320,Y(F) /S#180+200 

250 NEXT F 

400 LOCATE 3,20:PRINT "2=";2# 

410 LOCATE 3,21:PRINT "“A=";A 
20 LOCATE 3,22:FRINT "B=";E 
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Z=A+B/COSCT) 
Az 1 
B= 1 





3=A+B/COS(T) 
A= 3 
B= 1 






Z3=A+B/COS(T) 
A= 3 
B= 1 


Primitives et 
équations différentielles 


L'utilisateur propose une fonction y (x), l'ordinateur trace les gra- 
phes de y (x) et de deux primitives successives. Telle est l’idée première 
de ce programme, mais il permet également de tracer le graphe de la 
solution d'une équation différentielle du second ordre. 


Déroulement du programme 


La fonction y (x) est portée en ligne 1000. La fonction sera étudiée 
dans l'intervalle O, X. Dans cet intervalle, elle sera calculée pour 100 
valeurs de x (ligne 270). 

On peut éviter de recopier les lignes 10 à 40 qui ne constituent 
qu'un préambule. 

Pour lancer le programme lui-même, on fait RUN 200. 

La primitive de y est notée z (x), la primitive de z est notée t (x). 
L'utilisateur fixe X, z (O) et t (O). 


Pour chaque valeur de x, une valeur approchée de z (x) est calculée 
en ajoutant y (x).dx à la valeur précédente. L'ordinateur calcule de 
même t (x). Ces valeurs sont mémorisées. 


De plus l'ordinateur cherche si les valeurs de y, z et de t sont supé- 
rieures où inférieures aux valeurs précédentes. Ceci afin de déterminer 
l'échelle e en ordonnée de la représentation des graphes. 


Enfin il trace les courbes représentant y (x), z (x) et t (x). 
Critique 
Îl est bien évident que le pas d'intégration dx n'est pas infiniment 


petit et que donc les résultats seront approchés. Tel que le programme 
est écrit, dx est représenté par 6 pixels sur l'écran. 
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Vous pouvez le réduire à un pixel soit dx = X/600 en remplaçant 
100 par 600, 98 par 598, 6 par 1 en lignes 210, 270, 520 et 530. Vous 
pouvez également prendre un dx encore cent fois plus petit, mais en ne 
marquant qu'un point sur cent car un marquage plus serré serait inutile, 
et en ne mémorisant les valeurs de y qu'une fois sur cent pour ne pas 
surcharger la mémoire. 


Tout ceci améliore la précision du calcul mais augmente le temps 
nécessaire à la réalisation du programme... à vous de choisir. 


Détail du programme 


Lignes) 


10 présentation 
1000 fonction y (x) 


200 introduction des données 
280 situation pour x = 0 
320 boucle de calcul 


340 calcul de y (x) 

360 calcul de z (x) 

370 calcul de t (x) 

420 détermination des échelles 
470 représentation graphique 


Liste des variables 


x variable 

y (0,x) valeur de y (x) 

y (1,x) valeur dez (x) 

y (2x) valeur de t (x) 

f valeur maximale de x 

dx accroissement de x entre deux calculs 

k indice indiquant y ou zout 

n numéro du calcul (x = n.dx) 

p (k) valeur du maximum de y, zout 

m (k) valeur du minimum de y ouzout 

e(k) échelle en ordonnée dey,z out 

b (k) décalage en ordonnée pour y,zout 
afin de placer l'axe des abscisses. 
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Illustrations 


Figure 1 


t 


y = 4 avec z (0) = 0, t (0) 0 et X = 10 
la théorie prévoit z = 4.x et t=8.x 


Figure 2 


y = sin (x) avec z (0) = —1, t (0) = 0 et X = 10 
la théorie prévoit z = —cos (x) et t = —sin (x) 
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Figure 3 


y = sin (x) avec z (0) = 0, t (0) =0 et X = 10 
la théorie prévoit z = — cos (x) + 1 et t = — sin (x) + x 


Figure 4 


IF X > 10 AND X < 20 THEN Y = 1 ELSEY =0 
avec z (0) =0, t (0) =0 et X = 50 
Dans tous ces cas, l'ordinateur fournit des réponses satisfaisantes. 
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Figure 5 


t 


d?t{x) = 


Résolvons l'équation différentielle ax —t{x) 
x 





soit y = —t 

avec z (0) = 1, t (0) =0 et X = 10 

On s'attend à observer t (x) = sin (x) mais on note une divergence 
dans la réponse. C'est que le dx égal à X/100 soit 0,1 n'est pas très petit 
devant la période de sin (x) qui est égale à 2.PI. Le calcul est donc 
approché et devient très mauvais si l'on prend par exemple X — 100. 


Figure 6 


VSD. 





2 





Résolvons l'équation différentielle = sin (t (x)) 

soit y = sin (t) avec z (0) = 1, t (0) =1 et X = 10. 

On peut s'interroger sur le début de la courbe. Essayons alors 
X = 1 : le défaut a approximativement la même largeur, il est donc 


du au fait que notre procédé de calcul impose y (0) = O ce qui ne 
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correspond pas forcément à la solution de l'équation différentielle. 
Il ÿ a donc au départ de la courbe régime transitoire entre cette solu- 
tion forcée et la solution de l'équation différentielle. 


Figure 7 
LT 
Æ 
y 
d’?t 
Résolvons l'équation re = Xx.sin (t (x)) 
dx? 


soit y = x.sin (t) avec z (0) = 1, t (0) = 1 et X = 10 


Quel plaisir de pouvoir observer aussi rapidement la forme de cette 
équation différentielle ! 


Bien entendu toutes ces familles de courbes peuvent être utilisées 
comme exemples d'une fonction et de ses deux premières dérivées. Ceci 
permet de bien réaliser que la dérivée est liée à la pente de la fonction. 


Programme 


10 CLS 

29 LOCATE 14,6:FRINT"primitives" 

20 LOCATE 6,8:FRINT"et equations differentielles" 
49 LOCATE 1,22:PRINT'"'pour faire avancer le prograir: 
me, taper "CHR#(243) 

590 LOCATE 59,25:FRINT CHR$# (245) 

69 q=INKEY(1):1F q<>0 THEN 60 

70 CLS 

80 PRINT:PRINT"x varie de 0 a X (X>0)" 

90 PRINT:FRINT"z est la primitive de y" 

100 PRINT:PFRINT"t est la primitive de z " 

110 PRINT:FRINT'"'donner la fonction y(x,z,t)" 

120 PRINT:FRINT'"'en ligne 1000" 

130 PRINT:PRINT'"'puis faire RUN 200" 
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300 
910 
520 
530 
540 
5930 
569 
9570 
580 


LIST 1000 

REM 

CLS 

DIM y(2,100) 

FRINT:INPUT"'X ="; 

FRINT: INPUT'"valeur de z pour x=0 "“3;3y(1,0) 
PRINT: INPUT'"valeur de t pour x=0 ";y(2,0) 
PRINT:PRINT'"calcul en cours" 

REM 

dx=f/100 

x=0O: GOSUR 1000:y(0,0)=y 

FOR k=0 TO 2 
MHI=MIN(O,yY(k,0)):p(k)=MAX(0,y(Ck,0)) 
NEXT k 

FOR n=1i TO 99 

x=n#dx 

GOSUB 1000 

y (O,n)=y:p(0)=MAX(p(0) ,y):m(0)=MIN(M(O) ,y) 
y(i,n)=y(i,n-1)+dx#y(O,n-1) 

z=y(i,n) :p(1)=MAX(p(1),z2):m(1)=MIN(M(1) ,z) 
y(2,n)=y(2,n-1)+dx#y(1,n-1) 

t=y(2,n) :p(2)=MAX(p(2) ,t):m(2)=MIN m2) ,t) 
NEXT n 

REM 

FOR k=0 TO 2 

e(kK)=120/(p(k)-m(k)) 

b(k)=e(k)#m(k) 

NEXT k 

CLS 

LOCATE 1,4:FRINT"t" 

LOCATE 1,12:FRINT"z2" 

LOCATE 1,19:PRINT"y" 

FOR k=O TO 2 

MOVE 20,10+1350#k-b (k) : DRAWR 600,0 

FOR n=0 TO 98 

MOVE 20+6#%n,10+130#k-b (k)+e(kK)#y(Ck,n) 
DRAUR 6,e(k)#(y(k,n+1)-y(k,n)) 

NEXT n 

NEXT k 

END 

REM 


1100 RETURN 
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Equations différentielles 
du second ordre, complètes, linéaires 
et à coefficients constants 


Ce programme permet l'étude des équations linéaires, c'est-à-dire 
celles où d?x/dt?, dx/dt et x sont au premier degré, et dont les coef- 
ficients sont constants. Ces équations sont les plus utiles et les plus 
fréquentes dans les applications pratiques. 

Les équations étudiées sont de la forme : 

d? x dx 


de + Bon Cx + D = EsinWt 





A 


Le programme permet l'étude des oscillations libres de la forme : 
d? x 
dt? 





A + Cx=0 


Les oscillations forcées de la forme : 


d? 
A 7 + Cx= EsinWt 
dt2 


se : ; dx 
On peut aussi introduire une force retardatrice de la forme B re 


ou un terme de frottement sec. 


Nous rappelons qu'un «frottement sec » se caractérise par une force 
constante de signe opposé à la vitesse. 


Nous proposons à l'utilisateur un grand nombre d'exemples, il lui 
sera évidemment possible de changer les données pour traiter ses exem- 
ples particuliers. 


La résolution de l'équation différentielle est numérique, on peut 
donc reconstituer les régimes transitoires. 
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Liste des variables 


XO 
vo 


T 


ZOMOOX»> 


CX 


Dér 


abscisse initiale 
vitesse initiale 
coefficient du terme d?x/dt? 
coefficient du terme dx/dt 
coefficient du terme x 
introduction d’un terme permettant de simuler un frottement sec 
coefficient du terme sinWt 
pas du calcul 
durée de l'observation du phénomène 
et CT coefficients multiplicatifs de l'abscisse X et du temps T 


oulement du programme 


Lignes) 
10à 230 titre et entrées des données 


230 - 240 on trace les axes de coordonnées 
280 à 310 résolution numérique de l'équation différentielle et on 


trace la courbe sur l'écran 


340 si le coefficient CX choisi par l'utilisateur est trop grand, 
l'ordinateur avertit en inscrivant la phrase "CX trop 
" 
grand 


500 à 600 impression sur l'écran des différents paramètres 


10 
20 
30 
50 
60 
70 
80 
90 
100 
110 
120 
130 
140 
150 
160 
170 
180 
190 
Z00 
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Programme 


REM EGUATIONS DIFERENTIELLES 
REM 4462626 6260 2 DD DE DE DE DE DE DE DE DE DE DE DE DE DE 


REM 
REM EQUATIONS DE LA FORME: 
REM A#D2X/DT2+B#DX/DT+C#xX+D=EXxSIN (WT) 
CLS 
INFUT "VALEUR DE XO0="3, XO 
INPUT "VITESSE INITIALE VO="; VO 
INPUT "VALEUR DE A=";A 
INFUT "VALEUR DE E=";B 
INPUT "VALEUR DE C="3;C 
INFUT "VALEUR DE D=";D 
INFUT "VALEUR DE E=";E 
IF E:0 THEN INPUT "VALEUR DE W=";W 
IF E=0 THEN W=0 
INPUT "PAS DU CALCUL DT=";DT 
INFUT "DUREE DE L'ORSERVATION N='"3;N 
INFUT "COEFFICIENT DE X:CX=":;CX 
INFUT "COEFFICIENT DE T:CT=";CT 


220 MODE 2 

230 PLOT 0,0:DRAUWR 0,400 

240 PLOT 0,200:DRAWR 640,0 

250 LOCATE 2,2:PRINT "Xx" 

260 LOCATE 78,12:PRINT "T" 

270 X=XO0: TX=VO: M=ABS (XO) 

275 LOCATE 30,1:PRINT "XM="3;M:LOCATE 30,2: PRINT "T 
LD EE OO" 

280 FOR T=0 TO N STEP DT 

290 ATX=(EXSIN(WET)—-C#X—-EH#TX-—-D#SGN(TX))/A 
300 ACTX=ATX#DT: TX=TX+ACTX 

310 X=X+TX#DT 

3520 IF ABS(X):M THEN M=ABS(X):LOCATE 30,1: PRINT "X 
M="; M: LOCATE 30,2:PRINT "T="MST 

230 IF T#CT>=600 THEN GOTO 500 

340 IF X#CX>=180 THEN LOCATE 2,20:PRINT "CX TROP G 
RAND'":M=0: FOR Z=0 TO SOO:NEXT Z:GOT0 190 
350 FLOT T#CT,X#CX+200 

3560 NEXT T 

500 LOCATE 2,1S:FRINT “XO="3XO 

519 LOCATE 2,16:PRINT "VO="3; VO 

529 LOCATE 2,18:PRINT “A=";A 

530 LOCATE 2,19:FRINT "B="3:;B 

540 LOCATE 2,20:PRINT "C="3;C 

550 LOCATE 2,21:PRINT "D=";D 

560 LOCATE 2,22:PRINT "E=";E 

565 IF E<:0 THEN LOCATE 2,23:FRINT "W="3;W 
570 LOCATE 40,18:PRINT "DT="; DT 

580 LOCATE 40,19:PRINT UT=EUST 

590 LOCATE 40,20:PRINT "CX=";CX 

600 LOCATE 40,21:FRINT "CT=";CT 


XN= 6.944697326 
x 1: 0.33 





DT= 6.81 
T= 5.883993391 
QX= 158 


CI= 188 
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Méthode de Runge-Kutta 


Résolution d'équations différentielles 
par la méthode de Runge-Kutta 


Introduction 


Précisons tout de suite qu'il s’agit d’une résolution numérique, ce 
qui se traduira sur l'écran par une série de nombres. Si, par exemple, 
l'équation donnée, contient la variable X et la fonction Y (et ses déri- 
vés), on obtiendra des valeurs numériques de Ÿ en face des valeurs 
numériques de X. Dans toute équation différentielle, apparaissent des 
conditions initiales XO et YO. On donne alors à la variable un incrément 
H et l'on calcule une moyenne pondérée d’un nombre M d'évaluations, 
par la formule : 

Yi+1=Yi+{(K1+2K2+2K3 +K4)/6 


Les coefficients K1, K2, K3, K4 étant déterminés par : 


*K1=H.F (X;, Yi) 
*K2=H.F(X+1/2H,Yi+1/2K1) 
*K3=H.F(X+1/2H, Yi + 1/2 K2) 
*K4=H.F(X + H, Yi+K3) 


But du programme 
Cette méthode s'applique : 


1. À un système de deux équations du premier ordre de la forme : 
Y'=F(X, Y,2) 
Z'=G(X,Y,2) 


2. Aux équations du second ordre de la forme : 
Y"=F(X,Y, Y') 


En effet, il suffit de poser Y’ = Z, d'où Y” = Z' pour être ramené 
au premier cas avec : 


Y'=2Z 
Z'=F(X, Y,2) 
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Déroulement du programme 


Ligne (5) 
10à 165 Mise en page et entrée des données 


170 à 340 On calcule les coefficients K et on élabore les solutions. 
1000 et 1010 Le lecteur doit y inscrire les équations à traiter 


Dans un exemple, l'équation différentielle de départ est : 
Y"+4Y' +100 Y —100 X=0 


Si l'on pose Y’ = Z et Y” = Z', elle se transforme en : 
Z+4Z+100Y—100X=0 
ou Z'=—4Z—100 Y + 100 X, avec Y' = Z. 


C'est pourquoi, la ligne 1000 contiendra Y1 = Z 
la ligne 1010 ù Z1 =—42Z—100 Y + 100 X 


Après le RUN, l'ordinateur demande les valeurs initiales. 
Par exemple : X0 = 0.1 


Yo = 0.1 
Z0 =0 
H =0.05 
M =9 


qui permet d'obtenir 10 valeurs différentes de X et Y. 


X9= 6.1 
ÿ8= 8.1 

= 8 

= 8.85 

= 9 
X= 0.15 LVX= 6.181979167 ,2X- 8.114479167 
%= ERA Vr Er SSns" 7X= = 774169233 
X- 0.3 ,UX- 9.jg1c4am2 ,4X- 1126937239 
X= 0.35, VX- 6.254908296 ,2X- i. F7 0729 
X= 0.4 ,fX= 0. EE 15 ,2X= 1.51531458 
X= 6.45 ,VX- 6.484094386 ,2X- i.50557032 
X= 8.5 ,4X- eo 4726733056 ,X= 1.38549684 

= 0.55" ,YX- 0.541564272 ,7X- i. 20168693 
X= 8.6 ,ŸX= 0.526712437 ,2X- 1.89643196 


Y1=Z 

Z1=-4#2-1081+108xX 
eady 

CALL & 
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Autre exemple Retaper la ligne 1000 : Y1 = EXP (X) — Y 


90 

100 
115 
120 
130 
140 
150 
160 
165 
170 
180 
190 
195 
200 
210 
220 
230 
240 
250 
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Taper en 1010 :Z1=0 


Faites un RUN. Entrer XO = 0, YO —0, Z0=0, 


H=0.5, M = 4, par exemple. 


On a résolu cette fois l'équation du premier ordre 


Y' = EXP (X) — Y. 


1.1756847 
3.686954 


Y1=EXP(X)-Y 
218 


Ready 
CALL &h088 


REM METHODE DE RUNGE-KUTTA 
REM 36462606 6 DE DE DE DEEE DE DE DD DE DE DE HE HN 


CLS 
REM CONDITIONS INITIALES 
REM 3062606060 DD DE DEEE DE DE DE DE DD 
INPUT "VALEUR DE XO0="3; XO 
INPUT "VALEUR DE YO="; YO 
INPUT "VALEUR DE Z0=";2720 
INPUT "VALEUR DE L'INCREMENT H='";H 
INPUT "NOMBRE D'EVALUATIONS M=";M 
MODE 2 

LOCATE 1,1:PRINT “XO="3; XO 
LOCATE 1,2:PRINT "YO="3YO 
LOCATE 1,3:PRINT "ZO="3;20 
LOCATE 1,4:PRINT "H="3;:H 
LOCATE 1,5:PRINT "M="3;M 
PRINT 

FOR N=O TO M 
X=XO0:Y=YO:Z=270: GOSUB 1000 
K1i=H#Y1:L1=H#71 

X=XO+H/2: Y=YO+K1/2:Z=2Z0+L1/2 
GOSUB 1000 

K2Z2=H#Y1:L2=H#71 
Y=YO+K2/2:Z=Z0+L2/2 

GOSUB 1000 

KS=H#Y1:L3=H#71 

X=XO+H: Y=YO+K3S: Z=ZO0+LS 


8,521254225 ,ZX- 8 
3 ,2X= 8 

2,1202071$ 25 8 

="6.85336328 ,2{- 9 


260 
270 
280 
290 
S06 
5350 
340 
350 
9300 
310 
1000 
1010 
1020 
1030 
1040 
1050 
1060 


GOSUB 1000 
K4=HXY 1: LA4=H#Z1 
YX=YO+ (K1+2#K2+2#K3+K4) /6 
ZX=20+(L1+2#L2+2#L3+L4) /6 
FRINT'X=" SX" M "VX = VXG "M "ZX"; ZX 
XO=X : VO= VX: ZO=ZX 
NEXT N 
GOTO 1050 
REM SYST.EQUAT. (LIGNES 1000,1010) 
REM 62666 DEDE DEEE DE DE DE DE HE DE DE DEEE DEEE 
Y1SEXP (X)—Y 

Z1=0 

VH="EXP (X)—Y 

Wt= LL Ô LL 

RETURN 

LOCATE 1,20:FRINT "VI="; VE 
LOCATE 1,21:PRINT “Zi=";W# 
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Spirale de Cornu 


Problème 
1 
Soit à construire une courbe dont la courbure FR en un point soit 


proportionnelle à la distance S de ce point à un point O pris pour 
origine des arcs (distance comptée sur la courbe : arc S). On dit que : 


L ñ V? 
X=A cos 2 dV 
0 





V 2 





T 
et va f sin dV 
0 
En effet, posons : 7 V?=2U 
ji à 
X = — cos UdU 
: TV 


A 
va — sin UdU 
. ñmV 


A 
Donc le rayon de courbure R = —— et l'arc dS = R ; dU = AdV ; 
S = A.V. 7 V 


T 
Enfin R=—S 
A2 


La courbe cherchée est donc définie par les intégrales de Fresnel : 


: TV 
X=A cos 
0 


2 V 2 


TV 
dV v=a | sin 
0 


dV 








par 


Nous rencontrons une difficulté : celle de définir des coordonnées 








des intégrales. Pour la tourner, dérivons : 
dX ñV? 
— = À cos 
dU 2 
dY an V? 
— = A sin 
dU 2 


Système de deux équations différentielles que la méthode de Runge- 
Kutta nous permet de traiter. Pour la concordance des variables avec le 
programme, nous remplacons V par X 








X par Y 

Y par Z 
dY nX? 
— = 2 cos 
dx 2 ÿ : 

après avoir posé À — 

dz nr X2 
——= 2 sin 
dX 2 


Après RUN, l'ordinateur questionne : X0 0 YO=0 


H=0.02 M=1650. 


10 
20 
29 
30 
40 
930 
60 
70 
80 
970 
100 
110 
115 
120 
1:30 
140 
150 
160 
170 
175 
180 
190 
195 
200 
210 


Ligne 1000 par Y1 =2 cos (PI * (Y ? 2)/2) 
Ligne 1010 par Z1 —=2 sin (PI * (X ? 2)/2) 
Ligne 2030 : on trace la courbe. 


REM METHODE DE RUNGE-EUTTA 
REM SFIRALE DE CORNU 
REM 36263626 De DE DE DE DE DE DE DE DE DE DE DE DE D DE JE 
CLS 
REM CONDITIONS INITIALES 
REM 262626262626 26 De DE DE DE DE DE DE DE DE DE DE 
INFUT "VALEUR DE XO="3; XO 
INFUT "VALEUR DE YO="; YO 
INFUT “VALEUR DE ZO0=";70 
INFUT "VALEUR DE L'INCREMENT H=":H 
INFUT "NOMBRE D'EVALUATIONS M=":;M 
S=0:DIM A(M):DIM E(M) 

MODE 2 

LOCATE 1,1:FRINT “XO=":;XO 
LOCATE 1,Z2:PRINT "YO="; YO 
LOCATE 1,3:PRINT "ZO="3;2Z0 
LOCATE 1,4:PRINT “H="3;H 
LOCATE 1,S:FRINT "M=";M 
FOR N=O TO M 
LOCATE 1,7:PRINT N 
X=XO:Y=YO: Z=Z0: GOSUB 1000 
K1=H#Y1:L1=H#Z1 
X=XO+H/2:Y=YO+K1/2:Z=Z0+L1/2 
GOSUB 1000 
KH2=H#Y1:L2=H#Z1 


2 


Zo=0 
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220 
230 
240 
250 
260 
270 
280 
290 
‘300 
310 
32 

3530 
340 
ZS50 
300 
510 


Y=YO+K2/72:72=70+L2/2 

GOSUER 1000 

<Z=H#Y1:LZ=H#71 

X=XO+H: Y=YO+KZ: Z=20+L3 

GOSUB 1000 

H4=H#Y1:L4=Hx71 
YX=YO+(K1+2#KE2+2#K737+F4) /6 
ZX=20+(L1+2#L2+2#L7+L4)/6 
A(N)=YX:EH(N)=ZX 

IF ABS(YX):S THEN S=ABS (YX) 

1F AES(ZX)>S THEN S=AES(ZX) 

XO=X: YO=YX:ZO0=ZX 

NEXT N 

GOTO 10:30 

REM SYST.EOUAT. (LIGNES 1000,1010) 
REM EEE DEEE DE DEEE DEHÉ HE DE DE DE DE DE DE DE DE DE DE DE JE 


1000 Y1=2*#COS(FI#(X*2) /2) 


1010 


Z1=2*SIN(FI#(X"2)/2) 


1015 RETURN 


10750 


V+="2xCOS(FI#(X"2)/2)" 


1040 W#&="2#SIN(FI#(X"2)/2)" 


1050 
1060 
2000 
2010 
2020 


LOCATE 1,20:PRINT "Y1=";V# 
LOCATE 1,21:PRINT "Zi=";W# 
FLOT 320,0:DRAWR 0,400 
FLOT 0,200:DRAWR 640,0 
FOR N=0 TO M 


2050 FLOT A(N)/S#190+320,B(N) /S#190+200 
NEXT N 


2040 
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Y1=2xCOS(PI#(Xt2)/2) 
Z1=2xSIN(PI#(Xt2)/2) 





Y1=2#COS(PI#(X#2)/2) 
Don es) 


CALL 4088 


Le voyage terre-lune 
selon Jules Vernes 


En 1865, le très sérieux «Journal des Débats » publie en feuilleton : 
«De la Terre à la Lune», puis «Autour de la Lune» de Jules Verne. 
L'auteur imaginait alors la trajectoire d’un «boulet-wagon» entre les 
deux astres tiré à partir d’un gigantesque canon situé en Floride. 


En 1968, la nuit de Noël, le monde stupéfait suit à la radio et à la 
télévision les évolutions de la cabine Appollo 8, pilotée par F. Borman, 
J. Lovell et W. Anders, en orbite circulaire autour de la Lune. 


Quelques mois plus tard, le 16 Juillet 1969 à 3 h 56 mn 20 s (heure 
de Paris), le LEM (baptisé Eagle) de la mission Appollo 11 déposait 
les astronautes N. Armstrong et E. Aldrin sur le sol lunaire ; leur cama- 
rade M. Collins attendant leur retour en orbite lunaire dans la cabine 
Appollo. 


Un siècle plus tard, le rêve fou du Président Barbicane, du capitaine 
Nicholl et du journaliste Michel Ardan était devenu réalité. 


Le problème des trois corps 


Dans «Autour de la Lune», le Président Barbicane donne l'équa- 
tion du mouvement du boulet. 


Pour cela, il considère la Terre et la Lune immobiles dans un repère 
galiléen et applique le théorème de l'énergie cinétique entre les deux 
astres. L'équation obtenue est évidemment approchée. 


Nous avons dans un précédent ouvrage (1) étudié les interactions 
gravitationnelles et leurs applications à l'étude du mouvement des satel- 
lites terrestres. 


(1) «Programmes de physique sur AMSTRAD» - Editions Eyrolles. 
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Rappelons que : 


— La loi de gravitation exprime que deux masses M1 et M2 situées à 
une distance d l'une de l’autre, s'attirent mutuellement avec une force 
d'intensité : 


M1M2 
d? 





où G est la constante de gravitation universelle (G = 6.67 X 10° !! SI). 


Prenons le centre de la Terre pour origine des axes de coordonnées. 





Imaginons un vaisseau spatial de masse M en mouvement entre la 


Terre et la Lune. Au point de coordonnées X, Y à l'instant t, le vaisseau 
est soumis aux forces : 
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— 
— FT due à l'attraction de la Terre, 
— FL due à l'attraction de la Lune. 


— 


— — 
Ecrivons F = FT+FL 


Les composantes de F sont donc : 


G.MT.M G.ML.M 
FX = ————— EX) + : (XL — X) 
(X? + Y?) 2 (IXL=XT EN 72 
G.MT.M G.ML.M 
VE ———— EN ——————————— :" (—Y) 
NET (XE= RIT + YA 


— 
Appliquons la relation fondamentale de la Dynamique F = MY. Ù 
Les composantes de l'accélération du vaisseau sont donc : 


GMT GML 
VAR EX) + ——— (XL — X) 
a ds (XL — x)? + V2) 72 
GMT GML 
NE Ne A (Y) 
DE du (XL x)? + v2)5/2 


Liste des variantes 


R1 rayon terrestre au cube 
R2 rayon lunaire au cube 
MT masse de la Terre 

ML masse de la Lune 

XL distance Terre-Lune 


En ligne 460 la constante K est obtenue par l'opération suivante : 
400/XL (400 points sur l'écran correspondent à la distance Terre-Lune 
XL). 


Déroulement du programme 


Nous avons montré comment résoudre numériquement les équations 
différentielles. 


Pour calculer la trajectoire, nous utilisons l'algorithme de Runge- 
Kutta, précédemment exposé (voir la Spirale de Cornu), qui nous per- 
mettra une exploitation plus rapide du phénomène. 


Ligne (s) 


10à 260 Titres et introduction des données 
Le pas de calcul DT sera pris égal à 5 ou 10 mn 
300 à 490 On trace la Terre et la Lune (effet de loupe sur les deux 
astres par rapport à l'échelle de la trajectoire) 
500 à 580 Calcul de l'accélération 
900 à 1230 Utilisation de l'algorithme 
1250 à 2070 On trace la trajectoire 


Remarque : 
En lignes 530 et 540, le vaisseau percute la Terre ou la Lune. 


Les calculs effectués en appliquant tout simplement la relation 
fondamentale de la dynamique sont approximatifs, mais l'ordinateur 
nous permet de simuler le phénomène très rapidement (environ 20 mn 
par trajectoire). 


10 REM VOYAGE TERRE-LUNE 

20 REM MERE DEEE DE DE DE DEN NE 

30 REM 

40 CLS 

90 REM SIMULATION DE TRAJECTOIRES 
60 REM 366 EDEN DEEE DE DE DE DE DE DE EME DEEE EX 
70 KEM INTRODUCTION DES COORD.DE DEFART 
80 INPUT "X="3:X: INPUT "Y="sY 

90 KEM VITESSE INITIALE 

100 INPUT "VX="3; VX: INPUT "UVVY="sUvYy 
110 REM INTODUCTION DES CONSTANTES 
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910 

1000 
1910 
1020 
1050 
1040 
19060 
1070 
1080 
1090 
1100 
1110 
1120 
11:50 
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G=6.67E-11:MT=S.98E+24: MLe7.2SE+22: T0 
INFUT "FAS DU CALCUL DT=";DT 

CLS 

MODE 2 

LOCATE 5,17:FRINT "X="3;X3; "EM" 

LOCATE 3,18:FRINT "Y=";s EM" 

LOCATE 2,19:FRINT "VX="3; VX; "EM/S" 
LOCATE 3,20:FRINT "VY="SsUVYs"EM/S" 
LOCATE S5,21:FRINT "DT="; DT; "MN" 

Ri=2.5 D:R2=5.244E+18 














X=X#1000:Y=Y#1000:VX=VX#1000: VY=VY#1000 
DT=DT#60 

DIM A(4):DIM H(4):DIM C(4) 

DIM D(4):DIM E(4):DIM F(4):DIM N(4) 
LOCATE 20,2: FRINT"'VOYAGE TERRE-LUNE" 
REM TRACE DE LA TERRE 

REM 26262666 DEEE DE DE DE DE JE EDEN 

FOR I=0 TO 2:60 STEF 5 

J=I1*F1/180 

X1=60+20*#C0S (J):Y1=2O0#SIN(J) +220 

FLOT X1,Y1 

NEXT I 

REM TKACE DE LA LUNE 

REM 66 DD DE DE DEEE DE DE DE ED DE DE DEN 

FOR I1=0 TO 60 STEF S 

J=1+F1/180 

X1=460+10#C0S (J):Y1=10#SIN(J)+220 

FLOT X1,Y1 

NEXT I 

KH=400/XL 

XO=60+F#X:YO=220-F%Y:FLOT X0,YO 

GOSUEH 1000 

GOSUB 2000 

GOTO 475 

REM ACCELERATION DU VAISSEAU 

REM 6466 DEEE DE ED DE DE DE EME ED DE DE DE DE DE DEN 
Yi=H(1):Y2=H(2) 

N1=(Y1#YI1+Y2#V2)"1.S:1F Ni£R1 THEN STOF 
Z=((ABS(XL-Y1))"2+Y2#Y2)"1.S51F NESRE THEN 


N(1)=H(3):N(Z)=H(4) 
NCE)=-Ga#xMT#Y1/N1+GxMLx(XL-Y1) /N2 
NC) =-G#xMT#Y2/N1-G#xML#YZ/NE 
RE TUKN 
REM ALGORITHME DE RUNGE-EUTTA 
REM 62e 6 De DEEE DE DE DE DE DE DE DE DE DE DE DE DEEE DE 
A(1)=X:A(2)=Y:A(CI)=VX:A(A4)=VY 
FOR F=i TO 4 
R(F)=A(F) 
NEXT F 
GOSUEH 520 
FOR F=1 TO 4 
C(F)=N(F)#DT:H(F)=A(F)+C(F)/2 
NEXT F 
GOSUR 520 
FOK F=i TO 4 
D(F)=DTXN(F):H(F)=A(F)+D(F)/2 
NEXT F 
GOSUER 320 


1140 
1150 
1160 
1170 
1180 
1190 
1200 
1210 
1220 
12 
1250 
1260 
2000 
2010 
2020 
2050 
2060 
RES" 
2079 





FOR F=1 TO 4 
E(F)=DT4N(F):B(F)=A(F)+E(F) 

NEXT F 

GOSUR 520 

FOR F=1 TO 4 

F(F)=DT#N(F) 
A(F)=A(F)+(C(F)+2x(D(F)+E(F))+F(F))/6 
NEXT F 

X=A(1):Y=A(2): VX=A(Z): VY=A(4) 

RETURN 

REM TRACE DE LA TRAJECTOIRE 

REM #66 DEEE DE DEEE DEEE DE ÉD DE E 
X1=60+FxX:Y1=220-FXY 

IF X1*600 THEN GOTO 2050 

FLOT X1,Y1 

T=T+DT 

LOCATE 3,22:FRINT "DUREE-VOYAGE=";T/%6003% "HEU 


RETURN 


VOYAGE TERRE-LUNE 


CC) PRE EN NA Et Le O 


X= 6088 KA 
Ÿ= 8 KA 

11.5 KVS 
6ÿ= 6 QVS 


DT= 19 MN 

DUREE-VOYACE= 16.6666667 HEURESS 
Break in 548 
Ready 
CALL & 


VOYAGE TERRE-LUNE 


CD-—— æ | 


X= 8998 XN 
=8 
UX= 93.94 XVS 
UY= 6.12 KK/S 
DT= 18 MN 
DUREE-VOYAGE= 75.1666667 HEURESS 
Break in 58% 


Ready 
CALL &6898 
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VOYAGE TERRE-LUNE 


X= 5988 KA 

Ÿ= 0 KA 

UX= 9.915 KWS 

UY= 8.12 KVS 

DI= 10 MN 

DUREE-VOYAGE= 92.6666667 HEURESS 
Break in 578 


Ready 
CALL &h068 


X= 8006 KA 

Ÿ= 8 XX 

UX= 9.984 XA/S 

UY= 6.12 KW/S 

DT= 18 MN 

DUREE-VOYAGE= 165.833333 HEURESS 


Break in 5% 
Ready 
CALL & 


Calcul d'une intégrale définie 
Méthode de Simpson 


Introduction 


Etant donnée une fonction y = f (x) et sa représentation graphique. 





Soit à calculer l'aire comprise entre la courbe et l’axe des x, limitée 
aux verticales d'abscisses x, et x, . 
On partage l'intervalle x, —x, en N = 2p parties. 
Xp — 
N 


X 
D 
S f(x) dx est donnée numériquement par : 
X 


0 


. On démontre que l'aire cherchée : 





Soit E — 


[Yo + Yo) + 4 (ya + ya + + ypa) + 2 (2 + Ya + + Yp2)] 


œ |m 


S = 
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Déroulement du programme 


Lignes) 


10à 70 entrées des données XO, xP, N 
la donnée y = f (x) est introduite en ligne 90. Nous avons 
pris pour exemple y = 4/(1 + X T2) 
120 calcul de E 
150 on fait apparaître les termes figurant dans le crochet 
180 c'est le terme E/3 qui est multiplié par le crochet 


Exécution : RUN. Entrer X0O = O et XP = 1 N = 40. Vous retrour- 


10 
20 
50 
40 
50 
60 
70 
80 
90 
100 
110 
120 
1:50 
140 
150 
160 
170 
180 
185 
190 
200 
210 
220 
250 
240 


verez le nombre PI. 


REM METHODE DE SIMFSON 
REM 22626 26060636 DE DE DE DE DEEE DE 
REM 
CLS 

INFUT "HORNE INFERIEURE XO="3; XO 
INFUT “BORNE SUPERIEURE XP="; XP 
INFUT "NOMBRE DE DIVISIONS N="3;N 
LOCATE 12,12:FRINT "CALCUL EN COURS" 
DEF FNY(X)=4/ (1+X72) 

At="4/ (1+X°2) 

I=N/2-1 

E= (XP-X0) /N 

F1=0:F2=0: X=XO 

FOR E=1 TO I 

X=X+E:F1=F1+FNY (X) 
X=X+E:F2=F2+FNY (X) 


NEXT K 

R=E/3% (FNY CXO) +FNY (XF) +4% (F1+FNY (X+E) )+2%F2) 
CLS 

LOCATE Z,2:FRINT “Y=";A# 

LOCATE 3,4:FRINT "BORNE INF. X0="; XO 

LOCATE 3,S5:FRINT "RORNE SUF.XF="S: XF 

LOCATE %,6:FRINT "NOMB.DIV. N="3;N 

LOCATE 3,7:PRINT "REFONSE R='";R 

STOP 
Y=4/(1+Xt2) Y=LOG(X)/X 
PORNE INF.X8- 8 BORNE INF.X8- 3 
RON NS BORNE SUP.XP= 5 
NOMB .DIU. 48 NOAB.DIU. Ni= 46 
REPONSE R= Fa 14159266 RÉPONSE R= 8.691670715 

Break in 248 Break in 249 

Ready Ready 

CALL &4008 CALL &4098 


Calcul d'une intégrale définie 
Méthode des trapèzes 


Rappel du calcul d'une intégrale 
définie par la méthode des trapèzes 


Y 
Y=#(X) 


O X0 x1 XP 


Etant donné une fonction Y = f (X) et sa représentation graphique, 
soit à calculer l'aire comprise entre la courbe et l’axe des X limitée aux 
verticales d'abscisses XO et XP. 


On partage l'intervalle en N parties. On démontre que : 
XP 


(XP— X0) F1 1 
Î FAX). dX = ——— ETS + — F(XP) + f (X1) +.+#0P- 1] 
” N 2 2 
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Déroulement du programme 


1 


10 


30 
40 
90 
60 
70 
75 
80 
85 
90 
100 
110 
120 
150 
140 
1350 
195 
160 
161 
165 
170 
180 
190 
200 
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Le programme permet de calculer l'intégrale par intervalle choisi. 


Lignes) 


10à75 introduction des données 
80-85 on définit la fonction 
90 à 160 calcul de l'intégrale 
65 à 200 l'ordinateur imprime sur l'écran les différents paramètres 
nécessaires au calcul, ainsi que le résultat de l'intégration 


REM METHODE DES TRAFEZES 
REM 26026 JE DE DEEE 


REM 
CLS 
INFUT "BORNE INFERIEURE XO0="3:X0Q 
INPUT "BORNE SUPFERIEURE XP="3 XP 


INPUT "NOMBRE D'INTERVALLES N='":N 
LOCATE 12,1Z:PRINT "CALCUL EN COURS" 
DEF FNY(X)=4/(1+X"2) 
A$="4/(1+X"2)" 
E=(XP-XO) /N 

X=XO:R=O 


R= (FNY (XO) +FNY (XF) )/2 

FOR I=1 TO N-1 

X=X+E 

R=R+FNY CX) 

NEXT 1 

CLS 

R=R%E 

CLS 

LOCATE 3,2:FRINT A* 

LOCATE 3,4:FRINT "BORNE INF. XO="3 X0 

LOCATE 3,5:FRINT "BORNE SUP, XF="3; XF 

LOCATE 3,6:PRINT "NB.D'INTERV. N="3;N 

LOCATE 3,7:PRINT "REFONSE R=":R 
4/(1+Xte) LOGCX)/X 

RNE INF.XH- 9 Bo = 

BORNE SUP .XP= 1 BORNE S AU: Ve :XP= È 
NB.D'’INTERV. N= 48 NB.D: THEN, N- 46 
REPONSE R=- 3.14148849 REPONSE R- 6.6916679e2 

Ready Readu 

CALL &h0P8 CALL &h088 


Recherche des zéros d'un polynôme 
Méthode de Bairstow 


Soit p (x) =a, xn +a,., xn1 +...., +40,48n + 0 
un polynôme de degré n dont les coefficients 4,, 4,.,, … a, sont 
réels. 


Avec cette restriction, les racines de ce polynôme sont, soit réelles, 
soit complexes conjuguées. On peut donc l'écrire sous la forme : 


plx}=a, (x? +p;x+q,}1.(x2+p,x+q)72.... (x —x,) (x —x;) 
2m +m+.)+{(r+..….+i)=n 
Mettons ce polynôme sous la forme : 
p (x) = {x2 + px + q) p1 (x) +R (bp, g) x +S (p, q). 
Ceci est toujours possible puisque l’on a : 
(en posant : p, (x) =b, x"? +b, x +... + b;) 
bh = än 


bn-1 = 8n-1 — Pbn 


bg = 8x — pbk:1 — 4bk+2 avec K=n—2,....,0 


En supposant b, et b, définis de la même manière : 
R (p,q) = b; — a; — pb: — qb; 

S (p,q) = bo + pb1 = ao — gb = (ao — pb; — gb) + pb; 
En définissant des formules de récurrence analogues : 
Gb; 

Cn-1 = bn-1 — PCn 


C2 = b1 — pcs — qCa 
C1 = — PC — qC3 
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On démontre dans ces conditions : 
ôR ôR os ôS 


R,=—=—c Ro=—=—0c3 S)=—=c)—-pca Sa=—=—-0c — 
D 25 2 q 2q 3 D es 1 — PC2 q à 2 PC3 
La méthode de résolution est alors la suivante : 


— par une méthode itérative, chercher les valeurs de p et de qg qui annu- 
lent R (b, qg)et S (o, q). p et q seront ainsi les coefficients du trinôme 
que l'on peut mettre en facteur et dont les racines sont aussi racines de 
p (x). 


Cette méthode itérative est déduite de celle de Newton et elle s’ap- 
plique ici à un système de 2 équations à 2 variables*. 


On démontre que, après un certain nombre d'itérations, p, et q, 
tendent vers les valeurs de p et q cherchées : 


R 
Pa Pni SR (Pn-1, Qn-1) — RqS (Pn-1, Gn-1) 
RoSq — Sp Rq 
1 


Qn =Qn1 = = © SpR (Pn-1, Gn-1) + RpS (Pn-1, An-1) 
RoSg — SpRq 


— Connaissant ainsi les valeurs de p et de q, il est possible, en divisant le 
polynôme initial par ce trinôme, d'en déduire un polynôme de degré 
plus faible. 


— L'opération est réitérée jusqu’à l'obtention d’un polynôme de degré 
1ou2. 


Il est certain que la précision se dégrade au cours du calcul. 
* 1! faut que R (p, g) et S (p, q) s'annulent. 


En développant en série : 


ôR 

R (b,q) = R (po, Qo) + (p— 2e a De go) +{qg— go) — 3q (po, 4o) + … 
CSS 2S 

S (,q) =S (bo, go) + (bp — Po) 3 (po, Go) + (qg — go) (po, Go) + … 


Ce qui peut s'écrire : 


R (p, q) R (po, go) FU AE TU, go) P—Po 


(] 
+ 


S (p,q) S (bo, Go) 2e Po, ao) À Un, go) q—Qo 


En appliquant la méthode de Newton à ce système à 2 dimensions, 


on pose : 
R (p,q)=0 
S(p,q9)=0 


D'autre part, avec : 
ôR 


— ; )=R 
ee go p 


EE )=S 
dp 0, Qo D 


ôR 


rm o) = Rq 


os 
7 Po, 4o) = Sq 
ûqg 


Pour simplifier l'écriture, on a : 


R (po, Qo) Ro Ro P — Po 


S (Po, Go) Sp Sa q— Go 


D'où, après inversion de la matrice : 
-1 


P Po Rh Ro R 
. %o Sp Sa ? Po, go 
On sait que : 


1 _ 
Ro Ro Sp — Ra 
SO Sal |-S BR 
Rh Sq — Rg $ 
Partant du couple p,.:, q;.,, On trouve une nouvelle approche : 


-1 
Pn Pn-1 Ro Ra R 


An An-1 S Sal” Pn-1, Qn-1 


Rappelons que la méthode de Newton conduit à écrire : 
f(x) = f (xo) + (x — xo) * F' (xo) 
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qui 


10 
20 
FO 
40 
50 
60 
70 
80 
90 
» 9) 
100 
ine 
110 
120 
JT 
130 
#+ 
140 
159 
160 
179 
180 
190 
200 
210 
220 
230 
du 
240 
2350 
260 
270 
280 
290 
300 
310 
320 
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Puis, les racines de f (x) = O0 amènent à écrire : 
O=F (xo) + (x — x0) * F’ (xo) 

D'où x=x, HR 
f (xo) 

Soit l'algorithme général : 
F (x) 

F'Uxn) 





Xn+1 Xn — 


converge vers une racine réelle de f (x), si elle existe. 


Programme 


REM #66 260602 HE HE HE DE DEEE DE DE 


REM +*# # 
REM # racines d'’UN POLYNOME +*# 
REM +*# * 


REM 446266 DE DEEE HEDE DEN DE DE DE DE HE DEN 

REM d'apres W.lLuther 
REM Algorithmes 

REM editions Sybex 
REM ### ps:precision de solution du systeme r(p 
=0,s(p,q)=0### 

REM #*## pn:precision de la methode de newtont{l 
arisation)*### 

ps =0.0001:pn=0.000000001:ni=70:nj=ni 

REM #*## n degre du polynome p(x)=a(n)x"n+a(n-1 
(n-1)+...+a(0)##+# 

REM ##%# p(x)=(x"2+px+qg)#pl(x)+r(p,q)x+s(p,q) * 


INPUT “ degre de polynome n=";n 

PRINT #8," degre de polynome n=";n 

REM ###introduction des coeficients#** 

DIM an) ,bi{n) ,c(n) 

PRINT'"coeficients an) ,a(n-i),....a(0):" 

FOR k=n TO O STEP-1 

READ a(k):PRINT" “;a(k), 

PRINT #8,"";k,a(k) 

NEXT 

REM ### introduction des valeurs initiales p,q 


trinome x72+px+q ### 
PRINT 
INPUT" p="3;p 


INPUT" qg=";q 

PRINT:PFRINT "racines du polynome“ 

PRINT #8:PRINT #8,"racines du polynome‘ 
PRINT “methode de bairstow"” 

PRINT #8,"methode de bairstow“ 

1F n>2 THEN 510 

1F n=1i THEN PRINT-a(0)/a(i) 


250 IF n=il THEN PRINT #8,-a(0)/a(1):END 
340 pli=a(1i)/a(2):qi=a(0)/a(2) 

3250 re =-pl/2:im=rex*re-qil:lF im<O THEN 410 
360 FRINT re+SOR(im) 

370 FRINT #8,re+SOR(im) 

3580 FRINT re-SQR(im) 

390 FRINT #8,re-SQR(im) 

400 GOTO 450 

419 im=SQR(-im) 

420 FRINT resz;"+j#%";im 

450 FRINT #8,res;"+j#"s;im 

440 FRINT #8,re;"-jx";im 

450 n=n-2:1F n=0 THEN END 

460 FOR k =0 TO n 

470 a(k)=b(k+2) 

480 NEXT 

490 GOTO 3510 

500 REM +*##methode de bairstow ##+% 

9510 ni=0 

520 bin)=a{(n):c(n)=a(n) 

S30 bin-1)=a(tn-1)-pxb(n) 

540 cin-1)=b(n-1)-p#xc{(n) 

550 FOR k=n-2 TO © STEP -1 

560 b(kH)=a(k) -p*xb(k+1i) —-qxb(k+2) 

370 ckH)=b(k) -pxc(k+1) —-q#e(k+2) 

580 NEXT 

990 c(1)=c(1)-b(1) 

600 REM *##% iteration #*#*# 

619 d=c(2)#c(2)-c(1)*c(3) 

620 pi=p+(c(2)#b(1)-c(5)#b(0))/d 

68:90 qi=qt(c(2)*xb(0)-c(1)#b(1))/d 

640 IF (AES(b(1))<ps) AND (ABS (b(0)+p#b(1)) <ps) THE 
N 390 

650 IF (ABS (p1i-p)<pn) AND (ARS (qi-q)<pn) THEN 680 
660 ni=ni+i:IF niïnj THEN 700 

679 p=pi:q=ql:GOTO 530 

680 PFRINT'"precision";ps:;"non obtenue, choisissez d 
‘autres valeurs initiales" 

690 END 

700 FRINT'"pas de convergences apres'";nj;"iteration 
s,augmentez ni.'":END 

710 REM ##%% an), atn-i1),..., a(O)### 

720 DATA 1,9,13,0,26,0,-1 
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degre de polynome n= à 


8 1 
7 Q 
& [a] 
S (a) 
4 Q 
5 [a] 
2 [a] 
1 Q 
Ô —i 


racines du polynome 
methode de bairstouw 

1 
—1 
—1.209714E-07 +j% 1.00000012 
—1.Z20914E-07 —-j*x 1.0000001Z 
0.7071068575 +j+4 0.707106978& 
0.707106875 —-j+ 0.70710498& 
—0.707106887 +j+ 0.70710678S 
—0.707106887 —-j*# 0.7071046785 


degre de polynome n= S 
324 
—351 

108 

& 
—& 

1 


Cho bn 


racines du polynome 
methode de bairstouw 
0.33514458 
0.327235611 
0.334597249 +j%x 7.77077E-03 
0.334597249 —-j# 7.77077E-03 
—-0.:49999126 


p=i:q=i 


degre de polynome n= S 


S 324 
4 —351 
3 108 
2 & 
1 —8 
Ô 1 


racines du polynome 

methode de bairstouw 
0.333227791 +j*x 7.94569E-03 
0.333227791 -j* 7.94569E-03 
0.3447760SE& 

0.322312802 

—0.Z50000313 

p=-0.66 ; qg=0.111 

les racines exactes sont : 
1/3:1/3;:1/3,1/3,-1/4 


degre de polynome n= 6 


& 1 

1] 11 
4 45 
3 84 
2 70 
1 21 
[e) 1 


racines du polynome 
methode de bairstouw 
—0.50297S504 
—1.29079022 
-5.81176E-02 
—2.2410378S 
—3.13613972 
—-3.77090324 
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degre de polynome n= 6 


CHAT 


H Cm Cm 
ul 


mm © 
œ 


racines du polynome 
methode de bairstow 


0.197484 
-0.194256 


—5.26951E-10 
-5.249S1E-10 
1.:7590E-10 
1.27S390E-10 


+j+ 1.586955 
—j+ 1.566955 
+j+ 3.743350 
—j*+ 3.743350 


Statistiques et calcul 
de probabilité 
Loi de Gauss 


Introduction 


Le statisticien doit d’abord définir avec précision l'ensemble de 
référence sur lequel vont porter les observations. Cet ensemble s'appelle : 
« L'ensemble statistique ou population». 


Pour le physicien, le terme «population» désigne un ensemble de 
faits expérimentaux susceptibles d'être reproduits dans des conditions 
absolument semblables. 


Tous les facteurs possibles de variation qui n’entrent pas dans la 
définition de la «population», c'est-à-dire non contrôlés pendant 
l'expérimentation doivent être attribués au hasard. 


En fait, le calcul statistique ne s'applique, la plupart du temps, qu'à 
des populations dites «normales » représentées par la loi de Gauss. 


Il s’agit d'une loi de probabilité d'une variable aléatoire x continue, 
variable qui est définie de — © à +, 


Sa densité de probabilité est de la forme : 
1 {— X— Xmoy)? 
Vanot 9 20? 
Cette courbe est caractérisée par : 
— La moyenne arithmétique Xmoy des n résultats individuels X. 





FX) = 


— La variance : (X— Xmoy}?  — 
= = À 


y= —— = AXx? 


n 


avec X = X — Xmoy. 
— L'écart type o=Vv(E=0) 
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On peut définir aussi : 





100 
— Le coefficient de variation k«-= 
X moy 
; E 
— L'écart moyen ou erreur type : Emoy = — 
Vn 


— L'écart type affecte un résultat isolé et caractérise la précision de ce 
même résultat isolé. 


— L'écart moyen affecte une moyenne et caractérise la précision de la 
méthode analytique utilisée. 


En Physique : 


— Le terme «précision» désigne le degré de reproductibilité fonction 
seulement des erreurs accidentelles. 


— L'erreur systématique caractérise l'exactitude ou degré de justesse de 
la méthode. 


Les erreurs systématiques sont indépendantes du hasard et se pro- 
duisent toujours dans le même sens. 


Les erreurs accidentelles, au contraire, ont lieu dans les deux sens 
avec une probabilité égale. 


Loi de Gauss généralisée 


Une variable aléatoire continue X suit une loi de Gauss généralisée 
si elle a pour densité de probabilité : 


— (X— Xmoy) 


2 o? 





1 
f(X) = - 
(X) Vas exp 
La représentation graphique de f (X) est symétrique par rapport à 
(X — Xmoy). 








Loi de Gauss centrée réduite 
On passe de : 
f (x) . SE EN Rat (12/2) 
x) = RD = F= 
Vo? di 2 o? V2 Ars 


. X— Xmo 
par le changement de variable : + = +. 
Oo 
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Remarque : 
A l'occasion de ce changement de variable, le terme : 
X — moy 


7 07 ) devient exp (—t2?/2) et le terme constant : 
o 


exp —( 





1 
TH tiré de f (X) est multiplié par o. 

T0 

En effet, la division par © lors du changement de variable oblige 
pour conserver les surfaces et en conséquence les probabilités, à multi- 
plier les ordonnées par 0. 


Les problèmes qui concernent la loi de Gauss généralisée peuvent 
donc être résolus à l’aide de la loi de Gauss centrée réduite. 


— Lorsque t croît en valeur absolue, le terme exp (—t?/2) devient 
rapidement asymptotique à l'axe des abscisses. 


— La fonction f (t) passe un maximum pour t = 0 maximum égal à : 


1 1 
f (0) = —— 0) =——= 0,40 
(0) NZ exp (0) Ur 
— La fonction f (t) passe par deux points d’inflexion d'abscisses res- 
pectives — 1 et + 1. 


Déroulement du programme 1 


Ligne(s) 
10à70 entrées des constantes : 


. nombre d'intervalles : N 

. nombre de mesures : M 

. pas du programme : P 

. limite inférieure des classes : K 

. facteur d'échelle (dans le sens vertical) : F 


120-130-140 entrées des données numériques concernant la série 
étudiée. Nous avons rentré une série de mesures calori- 
métriques destinées à déterminer une constante d'étalon- 
nage 

150 à 450 l'ordinateur trace automatiquement l'histogramme de la 
série étudiée et indique la valeur moyenne 

500 à 560 on fait apparaître les différents paramètres nécessaires 
au bon déroulement du programme 
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10 
29 
25 


REM HISTOGRAMME D'UNE SERIE DE MESURES 
REM #22 DE DE DE DE DEN HE DH DE HD DE DE DE DE DE HE HE DE DE 


CLS: MODE 2 


30 INPUT "NOMBRE DE MESURES M=";M 

40 INPUT "NOMBRE D'INTERVALLES N='";N 

90 INPUT "PAS DU CALCUL F=";P 

60 INFUT "LIMITE INFERIEURE DES CLASSES K=";K 
70 INFUT "FACTEUR D'ECHELLE F=";F 

80 I=0:CLS 

85 LI=K 

90 FOR D=i TO N 

95 G=0: X=0:H=0: V=0 

100 FOR X=1 TO M 

119 READ Y 

120 DATA 2.0453,2.042,2.05:,2.067,2.029,2.065 
130 DATA 2.059,2.048,2.054,2.039,2.050,2.055 
140 DATA 2.030,2.042,2.027,2.0553,2.040,2.038 
150 H=Y+H 

160 IF Y>=K AND Y<=K+P THEN G=6+1 

170 NEXT X 

180 I1=1+G 

190 C1=G*#350/ (M/F) 

200 IF C12350 THEN GOTO 1000 

230 PLOT 0,0:DRAWR 400,0 

240 PLOT 0,0:DRAWR 0,400 

245 LOCATE 2,1:PRINT "P":LOCATE 50,24:PRINT "X" 
250 PLOT D#(300/(N+2)) ,0: DRAWR 0,G*#:50/(M/F) 
260 DRAWR 300/(N+2) ,0:DRAWR O,-(G*#350/ (M/F)) 
270 K=K+P:RESTORE 110 

280 NEXT D 

290 LOCATE 4,1:PRINT CINT(M/F) 

300 XMOY=H/M 

510 FOR 1=i TO M 

520 READ Y 

330 V=C(CY-H/M)#(Y-H/M)) +4 

540 NEXT I 

350 J=(((N+1)#P-(K-XMOY) ) #3500/ (N+2) ) /P 

360 PLOT J,0:DRAWR 0,400 

370 LOCATE CINT(J/8) ,1:PRINT "XMOY="; XMOY 
380 E=SGR(V/M) 

390 FOR I=CINT(SO00/(N+2)) TO 300 

400 X=P#((C(N+2)/300)#1-1)+LI 

410 T=(X-XMOY)/E 

420 G1=(1/SGR(Z#xPI))#EXF(-T#T/2) 

450 G2=(1/E)#(1/SQR(2#F1))#EXP (—(CX-XMOY) # (X—-XMOY) / 
(2#E°2)) 

440 PLOT 1,61#600:PLOT 1,2*6G2 

450 NEXT I 

900 LOCATE 45,3:PRINT "NBK.D'INTERVALLES=";N 
910 LOCATE 45,4:PRINT "NBR.DE MESURES="; M 

920 LOCATE 45,5:PRINT "LIMITE INF.DES CLASSES=";LI 
930 LOCATE 45,6:PRINT "LIMITE SUP.DES CLASSES=";LI 
+N#P 

540 LOCATE 45,7:PRINT "VALEUR MOYENNE=" 3; XMOY 
5930 LOCATE 45,8:PRINT "ECART TYFE='"3;E 

999 LOCATE 45,9:PRINT "FACTEUR D'ECHELLE F='"3sF 
560 STOP 

1000 LOCATE 53,3:PRINT "F TROP GRAND" 


1010 RESTORE 110 
1020 FOR W=1 TO 1000: NEXT W 
1030 GOTO 60 
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XMOY= 2.84633333 








AMEL: 
RÉESESS 
ES] 


XMOY= 2.84633333 
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Programme 2 


Le programme 2 permet d'étudier la fonction RND de l'AMSTRAD. 


5 REM ETUDE DE LA FONCTION RND 
6 REM DE DEEE EEE EE ED DE DE DE DE DE DE 


7 CLS:LOCATE 1,1:PRINT “APPUYER SUR A" 


8 A=RND 

9 AF=INKEYS 

10 1F A$<>"A" THEN GOTO 8 

25 CLS:MODE 2 

30 INPUT “NOMBRE DE MESURES M=";M:DIM A(M) 

40 INPUT “NOMBRE D'INTERVALLES N='";N 

30 INPUT "PAS DU CALCUL P=";P 

60 INPUT "LIMITE INFERIEURE DES CLASSES K=";K 
70 INPUT "FACTEUR D'ECHELLE F=";F 
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Lizk 

I=O:CLS 
FOR 
NEXT X 
FOR Dsi 


100 
110 
150 
160 
170 
180 
190 
200 
250 
240 
245 
250 


Y=A(X) 
H=Y+H 


NEXT X 
1=1+G 


TON 


G=0: X=0:H=0: V=0 
FOR X=1 TO M 


Ci1=G#350/ (M/F) 
IF C1>350 THEN GOTO 1000 

PLOT 0,0:DRAWR 400,0 

PLOT 0,0:DRAWR 0,400 

LOCATE 2,1:PRINT "P":LOCATE 50,24:PRINT "X" 
PLOT D#(300/ (N+2)) ,0: DRAWR 0,G#350/ (M/F) 


X=1 TO M:A(X) =RND 


IF Y>=K AND Y<=K+P THEN G=G+1 


260 DRAWR 300/(N+2) ,0: DRAWR 0,-(G#350/ (M/F) ) 


270 
280 
285 
290 
300 
310 
320 
330 
540 
350 
360 
370 
380 
390 
400 
410 
420 
430 
(2#E72)) 
440 
450 
9500 
510 
520 
530 
+NxP 
540 LOCATE 
550 LOCATE 
555 LOCATE 
960 STOF 

1000 
1910 
1020 
1030 
2000 
2010 
2020 


Y=A(F) 


NEXT I 


LOCATE 


NEXT I 
LOCATE 
LOCATE 
LOCATE 
LOCATE 
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K=K+P:RESTORE 110 
NEXT D 

IF I1<M THEN GOTO 2000 
LOCATE 4,1:PRINT CINT(M/F) 
XMOY=H/M 

FOR I1=i TO M 


CINT (J/8) 


E=SQGR (V/M) 
FOR I=CINT(300/ (N+2)) 
X=P#(C(N+2) /300) #1-1)+LI 
T=(X-XMOY)/E 
G1=(1/SQGR(2#P1))#EXP (-T#T/2) 

G2=(1/E)#(1/SQR(2#P1))#EXP (—(X-XMOY) # (X—-XMOY) / 


45,3: PRINT 
45,4:PRINT 
45,5:PRINT 
45,6: PRINT 


45,7: PRINT 
45,8: PRINT 
45,9:PRINT 


LOCATE 35,35:PRINT 
RESTORE 110 

FOR W=1 TO 1000: NEXT W 
GOTO 60 

LOCATE 3,53:PRINT 
FOR W=1 TO 1000 :NEXT W 
GOTO 40 


Va (CY-H/M) # CY-H/M) ) +V 


Jæ=(C(N+1)#P-(K—XMOY) ) #500/ (N+2) ) /P 
PLOT J,0:DRAWR 0,400 


°1:PRINT "XMOY="; XMOY 


TO 300 


PLOT 1,G1#600:PLOT 1,2*62 


"NBR. D'INTERVALLES=";N 
"NBR.DE MESURES=";M 

"LIMITE INF.DES CLASSES=";LI 
"LIMITE SUP.DES CLASSES=";LI 


“VALEUR MOYENNE=" ; XMOY 
"ECART TYPE=";E 
"FACTEUR D'ECHELLE F=";F 


“"F TROP GRAND" 


“N OÙ F MAL CHOISI" 


XMOY= 0.5B3524814 





8 
. CLASSES= 1.2 
HOYENNE- , . + is 


met 
AE 


[ay 
1 
Ce 





P 28 XMY= 6.536008492 





Programme 3 


Le programme 3 permet de calculer la population par intervalle 
choisi par application de la loi de Gauss centrée réduite (on utilise la 
méthode des trapèzes). 


10 REM LOI DE GAUSS 

20 REM 63e 

30 REM 

40 REM CACUL DE LA FOFULATION 

90 REM FAR INTERVALLE CHOISI 

6Q CLS:MODE 2 

79 INPUT "NOMBRE DE MESURES M=";M 
80 INPUT "VALEUR MOYENNE XMOY="3; XMOY 
90 INPUT "ECART TYPE E=";E 
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100 
N 

110 
120 
1350 


X—XMOY) / (2#E72)) 
R=(XP-X0O) /N 
X=X0:S=0 


140 
150 
160 
170 
180 
190 
200 
210 
220 
250 
240 
250 
260 
270 
O 
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INPUT 


"NOMBRE D'INTERVALLES (INTEGRATION) N="'"; 


INPUT "VALEUR INF.DE L'INTERVALLE XO0="3;X0 


INPUT 


“VALEUR SUP.DE L'INTERVALLE XF="; XP 


DEF FNY{(X)=(1/E)#(1/SOR(2%F1))xEXF (—(X—-XMOY) #« 


S=(FNY (XO) +FNY (XP) )/2 
FOR 1=1 TO N-i 


X=X+R 


S=S+FNY (X) 
NEXT I 


S=S#R 
FRINT 
PRINT 
PRINT 
PRINT 
STOP 


"VALEUR DE L'INTEGRALLE S=":S 
"POPULATION CORRESPONDANTE P="3;SxmM 
“POFULATION ARRONDIE FA=";CINT(S#M) 


REM ETUDE D'UN NOUVEL INTERVALLE:FAIRE GOTO 11 


NOMBRE DE MESURES M-? 18 
VALEUR MOYENNE XMO0Y-? 2.8463 

ECART TYPE E=? .011 

NOMBRE D'INTERVALLES CINTEÉCRATION) N-? 28 
VALEUR INF.DE L’INTERVALLE X9-? 2.81 

VALEUR SUP.DE L’INTERUALLE XP- 2.03 


VALEUR DE L’INTECRALLE S= 7.49536E-82 

SOPRATION CORRESPONDANTE P= 1,33296422 

POPULATION ARROMNDIE PA= 1 

Break HS 268 

Ready 

COT0 118 

VALEUR INF.DE L'INTERVALLE X6-2 4.83 
VALEUR SUP.DE L' INTERVALLE XP-2 2.84 
VALEUR DE L'’INTECRALLE S= 6. ENT 
POPULATION CORRESPONDANTE P=- 3.85199971 
POPULA NTI AREONDIE PA= 4 


eak in 2 


Ready 


co, de 
INF.DE L'INTERVALLE X0-? 2.84 
VALEUR SUP.DE L'INTERVALLE XP-? 2.85 


VALEUR DE L'INTECRALLE S= 8.339693382 
POPULATION CORRESPONDANTE P= 6.11449088 
POPULATION ARRONDIE PA= 6 

Break in 268 


A 
EE 
INF.DE L’INTERVALLE XB-? 2 
VALEUR SUP.DE L' INTERVALLE XP-9 £ 88 


UR DE L'INTECRALLE S= 6.997742948 
POPULATICN RRO TE D ANTE P 17.9593731 
POPULA FE ARRONDIE PA- 


Break in 
Ready 
CALL & 


Système d'équations linéaires 


411 X1 + 412 X2 és 
821 X1 + 822 X2 nn ee 


411 412 :....-... &in X1 b; 
43 422::.:.:.... 82n X2 b: 
X = 
dn1 4n2 *-...... nn Xn bh 
soit A.x = b. 


A est une matrice carrée n X n, x et b 2 vecteurs à n composantes. 


Il est facile de voir (figure ci-dessous) que, quand deux vecteurs sont 
orthogonaux, alors leur produit scalaire est nul : 
XX + Y, Y2 = 0 (= XY — YX). 
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Y2=-X 





Ceci se généralisant pour n dimensions, la relation (1) peut s’écrire : 


— 1.b; +811 x: +812 %2 + ns... + 81n Xn = 0 


— 1h +81 X1 + 8h2 X2 + ss. + 8hnn Xn = 


Ce qui revient à dire que le vecteur (— 1,x,,x:,....,x,) doit être 
orthogonal à tous les autres : 
(D, 811,812; 8inbss (b,,an1,:...8nn). 


On considère alors la matrice (n + 1) .{(n + 1) : 


db 63500 bn A (0) 
bD;,8@11,812, RS din A (1) 
bn, dn1,8n2,--...: »8nn A (n) 
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et l’on utilise l'algorithme de Schmidt qui permet de rendre la première 
ligne A (0) orthogonale à toutes les autres. Dans ce cas évidemment, la 
première ligne sera devenue : 


[1x1 x2...... Xn] (solution recherchée) 
sous la forme : 


[400 4801 ::.-.: 8on ] 


et donc : 


800 800 800 
Ce procédé consiste à normaliser la dernière ligne A (n) : 


A (n) 
[A (n)| 





—+ A (n) 


puis à remplacer A (n —1) : 
Aln—1) 


Aln—1)—[A(n).A(n—1)]e An + Ant) > —— 
n [A (n n ] n n AW 1) 


etc. 

Cet algorithme est répété jusqu'à ce que A (0) soit suffisamment 
orthogonal à tous les autres vecteurs. 

On détermine ensuite les résidus : 

r=b—A.x 

Posant A.x, = r (solution calculée de façon identique), il est cer- 

tain que x, + x est meilleure solution que x, , puisque : 
A.x; = b— Ax Alx+x;)=b 


Le programme 


Les données sont introduites dans le sous-programme en lignes 
610-700. Pour plus de commodité, elles se trouvent dans les DATA des 
lignes 1020 à 1100. 


Le programme proprement dit commence en ligne 720 (11 = nombre 
d'itérations). Une orthonormalisation est ensuite effectuée (sous- 
programme 70-510) et l'erreur conséquente estimée. 


Une première solution est alors estimée (800-830), un calcul des 
résidus est fait (840-950), leur valeur est affichée (930) et éventuelle- 
ment réutilisée pour une nouvelle boucle (ligne 110). 


La précision P peut être choisie égale à 1E — 9. 
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10 
20 
50 
4û 
50 
ca 


Le) 


NHNRNNE 
> OU BR ot) 
> © © © 


N N 
NN 
[Re] 


N 


ol 


ol 
= © Ÿ 
SG TGS 


Où 


b 
2 © 


y 
© 


oi O1 0! Où O1 0! 
in 
e) 


Q 
© 
[e] 


INPUT "ordre du systeme :";:N 
INPUT "precision :":P 
DIM Zz(N:N)D-a€CN:N) ;b (CN); x CN) >; sm(N) 
GOSUB S30 
GOSUE 720 
REM orthonormalisation 
FOR k=0 TO N 
FOR i= Q TON 
atk:i)=z(k:i) 
NEXT 
atk:Q0)=sm(k):a(0:k)=sm(k) 
NEXT 
ni=0 
FÜR g=n TO © STEP -1 
kKz=Q:m=0 
FOR r=n TQ g STEP -i 
ps=kz:kz=0:l=r+mim=i 
FÜR i=0 TO n 


ss=a(g:i)-a(l-i)*ps 
a{g-i)=ss:kz=ss*#a(r, i)+kz 
NEXT ii 

NEXT r 


IF kz=O THEN PRINT"det a=Q ?":END 
ps=1i/SCR(kz) 

FOR i = O0 TÜün 

a{g-i)=a(g>i)*ps 

NEXT :i 

NEXT g 

ps=0 

FOR g=i TO n 

kKz=0Q 

FOR i =0 TO n 

kz=a(0:i)#a(g>i)+kz 

NEXT 

ps=ABS (kz)+ps 

NEXT 

kz=ps/n 

PRINT"“orthogonalite solution=";:kz 
ni=ni+i 

IF (kz>p) AND (ni<4) THEN 140 

IF kz < p THEN 45Q 
PRINT'"precision “;ps;"non obtenue" 
IF AHS(a(Q:0Q)) > 1E-10 THEN 450 
PRINT"systeme incompatible ?7":PRINT 
IF a(0:0)=0 THEN END 

s=-1/a(0:0) 

REM solution 

FOR k=1i TO n 


490 
500 
510 
520 
9530 
340 
990 
560 
570 
580 
590 
6aa 
610 
620 
630 
640 
650 
660 
670 
680 
690 
700 
710 
720 
730 


740, 


750 
760 
770 
780 
790 
&ao 
810 
820 
830 
840 
6850 
860 
870 
880 
890 
900 
910 
920 
930 
940 
950 
960 
970 


a(G:k)=a(0:k)xs 
NEXT 

RETURN 

REM initialisation 
FOR k=1i TO n 

FOR i=i TO n 

READ z{k>i) 


NEXT 

NEXT 

FOR k=1 TO n 
READ b{(k) 
NEXT 


PRINT"'systeme d'equations lineaires" 
FOR k=1i TO n 

PRINT 

FOR i=1 TO n 

PRINT z(k,i)5" "3 

NEXT 

PRINT" ... "s3b(k) 

NEXT 

PRINT 

RETURN 

REM iteration 

ii=0 

FOR k=1 TO n 

sm(k)=b(k) 

x (k)=0 

NEXT 

REM execution d‘un pas 
GOSUB 70 

PRINT "solution apres";ii; "boucles" 
FOR k=1 TO n 

Xx(k) =x(k)+a(0;k) 

PRINT x{k) 

NEXT 

REM calcul des residus 
PRINT"residus" 

s1=0 

FOR k=1 TO n 

s=b{k) 

FOR i=1 TO n 
s=s-z(k;i)#xç(Ci) 

NEXT 

sm(k)=s 

PRINT s 

si=sit+s*s 

NEXT 

PRINT"moyenne quadratique:";s1/n 
IF s1/n < 1E-12 THEN END 
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eo 


1000 
1010 
1020 
1030 
1040 
1050 
1060 
1070 
1080 
1090 
1100 
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iisii+i 
990 IF 


ii > 3 AND s1/n > 0.00001 THEN 1010 
GOTO 780 

PRINT"“systeme incompatible ?":END 

REM matrice a par ligne 


DATA 
DATA 
DATA 
DATA 
DATA 
DATA 


74; 80; 18>:-11;,;-4,-8 
14,:-69,:21:28,0;7 
bb:-72»-5 7514 
—12:66;-30:-23; 3; -3 
3»:8>-7;-4; 1:0 
4:-12,:4:4,:0;,1 


REM second membre b 


DATA 


0:0:0;1:0:0 


systeme d'equations lineaires 


—74 80 18 —-11 -4 -8 ... © 
14 —-69 21 28 [e) 7 .. © 
66 —-72 -5 7 1 4 ... © 

-12 éé —-30 -23 3 -3 .. 1 
3 8 —7 -4 1 [e] .. © 


4 —-12 4 4 o 1 
orthogonalite solution= 
orthogonalite solution= 
solution apres O boucles 
—39.9999428 

34.9999468 

154.999771 
—1033.99848 
—3210.99526 

4095.99397 
residus 

1.52588E-05 
—1.S2588E-05 

[e) 

1.9073SE-05 

[e) 

—2.86102E-06 

moyenne quadratique: 
orthogonalite solution= 
orthogonalite solution= 


… © 
6.21609E-08 
9.70128E-11 


1.37607E-10 
4.43951E-04 
8.48862E-12 


Solution apres 1 boucles 
—40.0004501 

35.0003988 

155.001754 

—1034.0117 
—3211.03634 

4096.04635 

res i dus 

O 
—7.62939E-06 

[e] 
—7.62939E-06 

Oo 

[e) 
moyenne quadratique: 1.94026E-11 
orthogonalite solution= 7.37497E-04 
orthogonalite solution= 4.54747E-12 
solution apres 2 boucles 
—40.0001449 

35.0001241 

155.000556 
—1034.00371 
—-3211.01153 

4096.0147 
residus 


00CcO00o 


oO 
moyenne quadratique: © 


systeme d‘equations lineaires 
3 2 oo. 22 


2 -3 o.. 22 

orthogonalite solution= 2.79397E-09 
orthogonalite solution= 5.82077E-11 
solution apres O boucles 

1.69230769 

8.46153847 

resi dus 
—7.45058E-09 

[e) 
moyenne quadratique:z 2.77556E-17 
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80 


systeme d‘equations lineaires 
3 2 .. 22 


2 -3 ss 7 
orthogonalite solution= 1.92085E-09 
orthogonalite solution= © 
solution apres O boucles 

4 

5 
residus 

7.4S5058E-09 
—5.72529E-09 
moyenne quadratique: 3.46945E-17 


Changement de base 


Ce petit programme permet de connaître l'expression d'un nombre 
donné en base 10, en base 2 ou en base 6. |l est basé sur une suite de 
divisions. Dans le cas du passage en base 16, remarquez l'astuce qui 
consiste à utiliser le code ASCII du caractère pour appeler À, B, C... 


Voyons comment fonctionne le programme sur un exemple : 
conversion de 13 en base 2 : 


x=13 xx = INT (13/2) =6 
r = 3—2=1 
a$ = LA Es 
x — 6 
y$ = rar 

x= 6 xx = INT (6/2) =3 
r =66—-6-0 
a$ —_ 0” 
X — 3 
y$ == "01” 

x= 3 xx = INT(3/2)=1 
r = 3—2=1 
a$ _ ei 
x = 
y$ = "101” 

x= 1 xx = INT (1/2)=0 
r = 1 
a$ = LE 1 
x = 0 
y$ = 71101” 
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10 CLS 

20 PRINT:PRINT'"conversion d7/7un entier de base 10 e 
n base..." 

30 PRINT: INPUT"2 ou 16 "“;:b 

40 PRINT:INPUT'"nombre a convertir "3x 
50 IF x=0 THEN END 

SS IF bz=16 THEN GOTO 150 

60 REM 

70 xx=INT(x/2) 

80 r=x-2#xx 

90 a$=STRS(INT(r)) 

100 H=xXx 

110 y$=a$+ys 

120 IF x>1 THEN GOTO 70 

130 PRINT:PRINT Y$ 

140 REM 

145 END 

150 xx=INT(x/16) 

160 r=INT(x-16#xx) 

170 IF r<10 THEN a$=STR$(r) ELSE a$=CHR# (K+55) 
180 x=xx 

190 y$=a$+y$ 

200 IF x>1 THEN GOTO 150 

210 FRINT:PRINT Y# 

220 END 
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Logique 


Pour les bricoleurs en microinformatique, il y a un problème qui se 
pose extrêmement souvent et qui est celui du décodage d'adresse rela- 
tif à un circuit (périphérique) relié àune unitécentrale (microprocesseur). 
En effet, un microprocesseur 8 bits a un bus d’adresse constitué de 
16 fils ; cela correspond à environ 65000 possibilités ; chaque circuit 
relié à ce bus doit être sélectionné par une adresse (au moins) présente 
sur ces fils. 11 y a donc souvent une série de circuits logiques simples 
destinés à activer le périphérique (c'est-à-dire à lui appliquer «1» logi- 
que sur une entrée de sélection) lorsque l'adresse voulue est présente 
sur ces 16 fils. Le problème est donc : connaissant le circuit de décoda- 
ge, quelles sont les adresses possibles du boîtier ? 


En fait, un raisonnement souvent simple, permet de trouver la solu- 
tion ; n'est-il pas plus amusant de faire un petit programme ? l'algorith- 
me est simple ; il suffit d'examiner toutes les possibilités et de rechercher 
les bonnes solutions parmi elles. Enfin un programme qui fait faire à un 
engin idiot un travail idiot. 


Premier exemple (figure 1) 


Quelles sont les valeurs de A, ...... A, c'est-à-dire les adresses 
(sur 7 bits, donc 2 chiffres héxadécimaux) qui permettent d'obtenir 
L=tetM=1? 

On décompose le problème en 2 étapes : 


— chercher les valeurs de H, 1, J, K qui conviennent (lignes 40 à 140). 
Les différentes solutions sont rangées dans la matrice A. ZZ est leur 
nombre à la fin de cette étape. 


— pour ces ZZ solutions, on examine maintenant les valeurs de À, B, C, 
D, E, F, G qui y conduisent. Ce seront les adresses convenables. 


Remarque : 


— Porte NAND :L= NOT (A AND J AND |) 
— Porte OUEXCLUSIF :! =DXORE 
— Porte AND : AND 
— Porte OR : OR 
— Porte NOR : NOT (—- OR —- OR -) 
— "1" logique :—1 
"0" logique :0 


— La solution s'exprime en héxadécimal par : 


I=n 


D 
Il 
M 1 


2'.a; lici i =6) 


i 


0 


20 DIM A(30,30) 
40 FOR H=-1 TO © 

50 FOR J=-1 TO © 

60 FOR I=-1 TO 0 

70 FOR K=-1 TO © 

80 L=NOT(H AND J AND I) 

90 M=NOT (K) 

100 IF NOT(L AND M) THEN 110 

105 ZZ=ZZ+1:A(ZZ,1)=H 

106 A(ZZ,2)=J:A(ZZ,S)=1:A(ZZ,4)=K 
110 NEXT K 

120 NEXT 1 

130 NEXT J 

140 NEXT H 

150 REM SECONDE ETAPE 

160 FOR yy=1 TO 
170 FOR A=-1 TO 
180 FOR B=-1 TO 
190 FOR C=-1 TO 

200 FOR D=-1 TO 

210 FOR E=-1 TO 
20 FOR F=-1 TO 

230 FOR G6=-1 TO 

240 h=NOT( a AND b AND c) :i=d XOR e 

250 j=NOT(# OR g):k=)j 

260 IF (h=alyy,1)) AND (i=a(yy,2)) AND ( j=a(yy,3) 
) AND € k=atyy,4)) THEN r=64*ABS(g)+32*A0ES (f)+16%A 

BS (e) +8*ABS (d) +4*ABS (c)+2#ABS (b)+ABS (a): PRINT HEX# 
r,2)5" Me 

270 NEXT g:NEXT f:NEXT e:NEXT diNEXT c:NEXT b:NEXT 
a 

280 NEXT yy 


SOoSceoSon 


6R 25 4B 73 335 53 6D 2D 4D 75 55 55 69 29 49 71 31 
S1 6E ZE 4E 76 36 56 6A ZA 4A 72 32 52 6C 2C 4C 
74 54 54 68 28 48 70 30 50 7B 3K 5B 63 23 43 7D SD 
SD 65 25 45 79 59 59 61 21 41 7E SE SE 66 26 46 
7A SA SA 62 22 42 7C SC SC 64 24 44 78 58 58 60 20 

40 6F 2F 4F 77 57 57 7F SF SF 67 27 47 
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Fig. 1. 


Second exemple 


Celui-ci est extrait de la notice d’un kit d'évaluation pour micro- 
processeur de la firme HEATHKIT (marque déposée). Les lignes d’adres- 
se À,: à À, sont reliées au décodeur d'adresse constitué par 3 portes 
NOR (IC — 2A, IC — 2B, IC — 2C) et une porte NAND à 8 entrées 
(IC — 3). Les lignes À,, A,, À, permettent de sélectionner 7 registres 
internes de IC — 1. Celui-ci est validé par E = O, donc par une sortie de 
IC — 1 à 0” (appelant K cette sortie, il faut donc avoir K = 1). 


Le programme évalue d’abord les sorties de IC — 2, appelées |, H, 
G ; K est aussitôt évalué et s'il vaut ”1”, la solution est affichée en 
héxadécimal. Une adresse seule convient ici : c'est 02F0 ; les 8 registres 
internes sont donc validés par les adresses 02F0 à 02F7. 
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aN9 "& ‘614 





z£l 
-€vv 


O€S Tv 
€-91 





G+ 


0€S ve 
€-91 
















pl 
SX2018 

nie JOVAY3INI 
JOONY pl si SS3400v 

Nas WINIVUL 

2N 21 af 01 

6 1! CE 0000 0000 

> 01 : À Ovlvtvév Sy vlv VAT NN 

40 6 ç+ V  Y z— 3N1 

FE 000 0000!sY018 32VIU31NI 0400 

à 1 { [0000 0000 VLVO YINIVYI 0000 

ae Ualacat vaSa°at EME TU 

2N RDA Ne ggevwa "HI 

2N » 1 
J0ONv € UMAIA LNOU3) 19 S PE 2 1 0 SS3YO0v 

j 2 da 6 1 à p 2 q € 1N3W93S 

ee: 

Nid 


10 FOR 
—1 TO 0 
20 FOR 
30 FOR 
40 FOR 
9590 FOR 
60 FOR 
70 FOR 
80 FOR 
90 FOR 
110 FOR 
120 FOR 
130 FOR 
140 FOR 
150 FOR 
160 I=N 
170 H=N 
180 G=N 
190 K=A 
ND I 
200 IF 
AS)+16* 
56*ABS ( 
6*ABS (C 
NT HEX$ 
210 NEX 
220 NEX 
230 NEX 
240 NEX 
250 NEX 
260 NEX 
270 NEX 
280 NEX 
290 NEX 
300 NEX 
310 NEX 
320 NEX 
339 NEX 
340 NEX 


RUN 


02F0 


AS= 


A4=-1 TO 
AS=-1 TO 
A6=-1 TO 
A7=-1 TO 
A8=-1 TO 
A9=-1 TO 
s-1 TO 0 
=-1i TO © 
C=-1 TO O 

D=-1 TO © 

E=-1 TO © 

F=-1 TO © 

VMA =-1 TO O0 

OT(E OR F OR VMA ) 

OT(B OR C OR D) 

OT(A3 OR AB OR A) 

4 AND AS AND A6 AND A7 AND A9 AND G AND H A 


CGCoOCcCoco 


K THEN R=ABS (AO) +2*ABS (A1) +4#ABS (A2) +8#ABS ( 
ABS (A4) +32*ABS (AS) +64#ABS (A6) +128*ABS (A7) +2 
A8) +512*ABS (A9) +1024%ABS (A) +2048*ABS (B) +409 
)+8192*ABS (D) +16384#ABS (E) +32768*ABS (F) : PRI 
(R,4) 

T VMA 

TF 
TE 
TD 
TC 
TB 
TA 
T A9 
T AB 
T A7 
T 46 
T AS 
T A4 
T A3 


Troisième exemple {figure 3) 


Quelles sont les adresses possibles de IC — 1 ? (Attention : À,, à 
A, ne sont pas employées). 
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SX2018 
39V3431NI 

SS3#00v 
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0000 0000 
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O8I|_[VWA 
Le] [e] 
0000 


avai Liu 


Hans) 0000 
0alaata vaSsa9aia 


SX9018 39VIY3INI 
VLVO YINIVYL 





oao!3y 


[ 
0/1 vivo 
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Séries de Fourier 


Un signal périodique, de période T, est décomposable en série de 
Fourier : 


f(t)=ao + Z a, cos nwt +b, sinnwt 


T 


’ 
0 = + f f(t) dt 
T 
0 
2 T 
8n — A f (t) cos nwt dt 
* 0 
T 


2 
bh =+f f (6) sin not dt 
0 


Si l’on connaît la période T, il est possible de calculer numérique- 
ment ces intégrales. 
| : T 2T 
Si f (t) est connue sur N points f (0), f = ,f (©) ose 


(N—I)T\ . 
f|————}, l'intégrale donnant a, devient : 
2 isNA ff T\ T 
an = — 2 {£) cos (no =) ST 
T i=0 N N N 


Pour accélérer les calculs, on peut poser : 


noT n2rT 27n 


89 


Nous allons utiliser cette méthode pour calculer les séries de Fourier 
de quelques signaux. Nous découvrirons au passage les limites d’utilisa- 
tion de cette méthode et la manière d'interpréter les résultats qu'elle 
fournit. 


Tous les programmes proposés sont calqués sur le même modèle, 
aussi ne commenterons-nous que le premier complètement. 


Les deux premiers concernent le calcul d’une série de Fourier : 
calcul numérique et représentation graphique d’un spectre. 


Le troisième montre l'influence sur le spectre d’un signal de la trans- 
mittance d'un circuit. 


Les quatrième et cinquième programmes démontrent l’utilisation 
qui peut être faite d’une série de Fourier : évaluer graphiquement la 
tension de sortie d’un filtre passe-bas du premier et du second ordre. 


Programme 1 


n points de calcul 

nh nombre d'harmoniques calculés 

j variable de boucle (varie de O à nh) 

8, coefficient du terme en cosinus correspondant à l’harmonique 
b, coefficient du terme en sinus correspondant à l’harmonique j 


Il permet d'afficher les valeurs numériques des coefficients a, et 
b,. Nous donnons à titre de comparaison les valeurs exactes de ces 
coefficients pour un signal carré. On peut également étudier le module 
de chaque terme, soit SQR (22? + b2). 


Pour un signal carré, on a : 


T 
f{t) = 1 ROÉDRRRSS 


T 
f(t) =—1 POUR t<T 
an = 0 quelque soit n 
4 
bh = — pour n impair, sinon 0. 
mn 


On constate que la précision augmente quand on prend un plus 
grand nombre de termes (n grand). 


90 


10 REM  #etkstesk sets teteote se ef ee ee he ee 


20 REM * # 
30 REM * serie de fourier * 

40 REM * * 

50 REM * programme {À + 

&Q REM * * 

70 REM  #3#k3#44e0k0#e ee ie eee RO 

80 KREM 

90 CLS 

100 PRINT:FRINT'entrez la fonction en ligne 110" 
110 LEF FNfCi)=SSN(SIN(Z*#FI#i/n)) 

120 INPUT " nombre de points de calculs"sn 
130 INPUT "harmoniques calcules :"“;nh 

140 FOR j=0 TO rh 

150 an=0:bn=0 

160 a=2*P1IxJ/n 

170 FOR i=0 TO n-1 

180 an=FNfCi)*CiS (ai) +an 

190 bn=FNfCi)kSIN (aki) +bn 

200 NEXT :i 

210 bn=£#bn/n 

220 an=#+an/n 

230 GOTO 250 


240 
230 
z60 


PRINT ‘"a";J;"=";an; "b"sJ;"=";:bn 
PRINT “nh.."s js". .—> "3 SGR (an £+bn7"2) 
NEXT j 


FNfCI)ESGN(SIN(Z*FI#i /1) ) 
nombre de points de calcul: 3 
harmoniques calcules : 10 


hs © > Q 

nh.s 1 ...—> .2671463 
vhs. 2 ...—> S.57056E-10 
mhes S > 0.41:2069776 
rh. 4 ...—> 4.3:173E-10 
nh., S > O0.2338523S552 
nh.. 6 ...—> 1.1976S4E-07 
rhos 7 e.—? O0.15235129741 
rh  ...—> 4.11590E-11 
nhos 7 > O0.102284547 
rh. 10 ,..—> £#.18279E-11 
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FNFCIi)=SGN(SINCEKPIKi /m)) 
nombre de points de caleuls 1##% 
harmoniques calcules : 10 
Ms Où eue [a] 

he. 1 ...—> 1.27298357 
mhes LE uso 3.14317E-10 
nbss So ous > O.4:53445914 
mh.. 4 ...—> 1.40575E-10 
he 3 > O,243345307 
mh.. 6 ...—> 1.10900E-07 
nhes 7 cu. Q.15009S188 
vhs 3 oo "JE-10 
nh.s 9 ...—à 

mh.. 10 ...—7 
valeurs exactes 

oh 1 ...-r 1.27323754 
..nh 3 ...-5 0.4244131%* 
.nh 5 ...—r O0.254447907 
th 7 ...—> 0.1#18715344 
..nh 7 ...—> 0.141471041 
oh 11 ...—> 0.11574708S 


Programme 2 


il est identique au précédent jusqu'à la ligne 260. A partir de là, 
l'amplitude des harmoniques est stockée dans une matrice Y (nh). La 
valeur maximum mx des coefficient est évaluée. Puis le spectre du signal 
est tracé tout en étant normalisé (la plus grande raie aura toujours 
400 points de hauteur). 


Nous proposons un exemple particulier et pour lequel théorique- 
ment le théorème de Fourier ne s'applique pas. En effet, sur la périodeT, 
nous avons construit (ligne 110) une fonction sinusoïdale de période 
T 
56 Comme le calcul ne peut mettre en évidence que les coefficients 

; 2 3 
correspondant aux fréquences — ,—,—,.. 
P jé TT 


nous donner une raie à a . Pourtant, on constate que le spectre (faux) 


.., il ne peut évidemment 


2,5 ; 
obtenu est à peu près centré sur TT et que finalement ces résultats 


sont exploitables. Nous reviendrons plus loin sur ce phénomène. 
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D'autre part, et seconde constatation indépendante de la précédente, 
que se passe-t-il si l’on cherche à évaluer les harmoniques «lointaines », 
pour lesquels le calcul ne devrait pas être fait ; en effet, il y a alors trop 
peu de points et l'équivalence : calcul de l'intégrale = calcul numérique 
n'est plus justifiée. La seconde figure montre le résultat : le spectre 


n 
obtenu est symétrique par rapport à FT Résultat remarquable ! En 


effet, sur une période T, nous prenons n valeurs de la fonction étudiée. 
Cela suppose implicitement que cette fonction est nulle le reste du 
temps. 


Conclusion : nous pensons effectuer le calcul approché des termes du 
développement, mais nous effectuons le calculexact de ces termes pour 
une fonction définie seulement en quelques points. Et le spectre obtenu 
est justement celui-là, qui ressemble dans sa moitié gauche à celui qui 
est cherché. 


Il y a ce que l'on appelle «recouvrement » ; c'est-à-dire que les raies 
; n ; : 
de fréquence supérieure à PT viennent se mélanger aux raies de fré- 


n 
quence inférieure. En conclusion, le résultat (le spectre de O à 2T } ne 


peut être exact que si le signal de départ n'a pas d’harmoniques de fré- 


n 
quence supérieure à 2T : c’est le théorème de Shannon. 


10 REM #64 AO AO UE 
ZO REM # 
30 REM 
40 KEM 


# 

# serie de fourier # 

# 
50 REM # programme Z # 

4: 

# 


# 


60 REM #4: 
79 KEM 
80 KEM 
90 CLS 
100 PRINT'entrez la fonction en ligne 110" 
110 DEF FNfCi)=SIN(E#PI#i#2.5/n) 

120 mx=0Q 

130 INPUT " nombre de points de calculs"sin 
140 INPUT "harmoniques calcules :":nh 

150 DIM yCmh) 

14Q FÜR j=0 TO nh 

170 an=Othn=o 

150 a=#PIl#J/n 

190 PRINT "j="; ) 


ee 24e ee ee AC AH EME 
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oi Gi ti 


o1 D 5 


FOR i=O TO n-i 

an=FNf Ci) Ca+i) +an 
bn=FNfCi)#SIN Caxi) +bn 
NEXT :i 

bn=z+tbn/n 

an=£#an/n 

ÿCj) IR C(an+tant+tbrn:#ttn) 

IF J=0 THEN 70 

IF y(j):y(j-1) THEN mx=y(j) 
NEXT j 

MODE 

FLOT O:Z:DURAW 440,2 

FÜR J=0 TO nh 

FLOT J#439/vbh: OS: DRAWR Q:yCJ)#3%%/mx 
NEXT J 





FNf{1)=SIN(DaPlIxi#z.S/n) 
nombre de pointe de calcul: 32 
harmoniques calcules : 10 


FNFCIi)=SIN(Z+FIxi#2.S/n) 


nombre de points de calcul: 32 


harmoniques calsules : 32 


Programme 3 


Il montre l'exploitation qui peut être faite du théorème de Fourier. 
Un signal (carré en l'occurence), dont le spectre est calculé (et n'est 


4 LA LA n e . . 
retenu que jusqu’à l’ordre 3 ligne 150 pour les raisons qui précèdent), 


traverse un filtre RC du premier ordre, dont la transmittance (lignes 
370 - 390) est tracée point par point. Le programme déduit le spectre 


du signal de sortie, (légèrement décalé pour des raisons de visibilité). 


10 

zO 

50 

40 

50 

COTE) 

70 

30 

90 

100 
110 
120 
130 
140 
150 
140 
179 
130 
1970 
200 
Z10 
220 
z'30 
40 
£50 
Z26Q 
270 
250 
270 
300 
310 
320 
330 


REM Hits te secte ce fe este he he Re 


REM *# # 
REM # serie de fourier * 
RENM * # 
REM * programme 3 # 
REM # 4e 
REM #4 te et M ER AH 
REM 

CLS 


PFRINT'erntrez la fonction en ligne 110" 
DEF FNf(i)=SGN(SIN(Z2#PI#i /n)) 


mx=(Q 

INFUT " mombre de points de calcul:"sn 
INPUT "rapport f/fce :"5k 

nh=n/2 

DIM y(nh) 


FÜR j=0 TO nh 

an=0: bn=0 

a=2*+PI+J/n 

FRINT "j="3) 

FOR i=09 TO n-1 
an=FNfCi)xC0S (ai) +an 
bn=FNfCi)#SIN (a+i) +bn 
NEXT ji 

bn=Z#bn/n 

an=Zz#an/n 
ÿY(j)=5SGR(an#xant+tbn#En) 
IF y(j)>mx THEN mx=y(j) 
NEXT j 

MODE Z 

PLOT OQ,2:DRAW 640,2 

FOR J=0 TO nh 
s=y(j)/SGR(1+(j#xk) T2) 
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340 
350 


340 


370 





NEXT J 
FOR J=0 TO nh 
PLOT J#600/nh: 379/SER (1+ Cj#k) 2) 


370 


FLAT (J+0.5)#€00/vnh; O3: DRAWR O:s#377/mx 


FLOT J#6400/nh; OS: DRAUWR 0:y (CJ)+3577/mx 


STEF nh/400 


NEXT 


FN CI) ESGNCSINECZ#FI#i /n)) 


nombre de points de 
5 1 


calcul: 44 


rapport f/f2 





FNfCI)=SGNCSINCE#PFI#i /n)) 
“ombre de points de calzsul: 44 


rapport f/fc : O©0.Z 


Programme 4 


Il permet, dans les mêmes conditions, de tracer la courbe donnant 
l'aspect du signal de sortie de ce filtre, en tenant compte de l’atténuation 
et du déphasage introduits par celui-ci sur chaque harmonique. 


10 REM #4 et A A HER 
ZOO REM *# 


30 KEM # serie de fourier # 
40 REM * # 
50 REM # FROGRKAMME 4 # 
4O REM * # 
70 REM #24 RH MO HER 


50 REM 

90 CLS 

100 FRINT'entrez la fonction en Ligne 110" 
110 DEF FNfCi)=SGN(SIN(E#PI#i /n)) 

120 REM -# la fonction doit etre impaire 
130 mx=0 

140 INPUT " nombre de points de calculri"z:n 
150 INPUT "rapport f/f2 5"5kK 

160 nhen/?2 

170 DIM y(nh) 

180 FÜR j=0 TO vh 

an=03 bn=Q 

a=£2#PI#gJ/n 

FRINT "j=";)j 

FOR i=Q0 TÜ n-i 

bnzFNfCi)#SIN (ai) +bn 

NEXT ii 

bn=2+bn/n 

y(j)=En 

NEXT j 

MODE 2 

FPLOT O:ZO0:IDRAW 640; #00 

oO FÜR t=0 TO 600 

© FÜR J=0 TO vnh 

S=y (j) SIN CEKPIX j4#t/GOO-ATNCj#E) ) /SER (14 Cj#E) TE) 
O ss=ss+s 

NEXT j 

» FLOT t:200+100%55 

540 Sss=0 

D NEXT t 





0) 


D 4 4 
Di Hi 
Le] 


Di 0j 
on à 
c © 


ü 
D 
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nombre de points de calcul: 64 
rapport f/fc : 1 


nombre de points de calcul: 64 


rapport f/fc : O.1 





Programme 5 


Il réalise la même opération pour un filtre du second ordre, pour 
lequel le rapport F/F, (fondamental sur fréquence naturelle) et le coef- 
ficient d'amortissement m sont donnés. 


Si, avec la méthode utilisée dans les programmes précédents, et avec 
les précautions d'emploi énumérées, il est possible d'évaluer une série de 
Fourier, il n'en reste pas moins que le temps de calcul devient vite 
prohibitif dès que n est assez grand, disons supérieur à 64. En effet, 
doubler n permet d'évaluer à 2 fois plus de raies, au prix d’un plus grand 
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nombre de multiplications pour l'évaluation de chaque intégrale. En 
fait, à n donné, il faut 7 multiplications par raie du spectre, d’où environ 
un temps de calcul proportionnel à n2. 


Deux mathématiciens, Cooley et Tukey (1965) ont réussi à optimi- 
ser ces calculs et ont mis au point un algorithme, applicable lorsque n 
est une puissance de 2, dont le temps de calcul est proportionnel à 
nlogn, en se servant de symétries qui interviennent pour ce choix de n. 
Cette méthode sera développée dans le chapitre suivant. 


10 REM Hit 
“O0 REM 
30 REM 
40 REM 


# 

+ serie de fourier 

# 
50 REM * FROGRAMME © 

#: 

# 


4 
x: 
# 
40 REM : 
70 REM 
50 REM 
90 LS 
100 FRINT'"'entrez la fonction en ligne 110" 
110 DEF FNfCi)=ESSNCSIN(E#FI#i /n)) 
130 REM —* La fonction doit etre impaire 


+ 
D ER ER LP LEE ED ER D EE 5 ED 9 5 25 LE 


130 mx=0 

140 INFUT " nombre de points de calculi"s;n 
150 INFUT "rapport f/fc s"3sEk 

1640 INPUT "coeff. d amortissement:";m 
179 nhen/# 

130 DIM y (nb) 

190 DIM t(nh):LIM p(nh) 

200 FÜR j=0 TO vh 

10 rat-(Cj#k) Es ii zZémej#k 

LEO t(j)=1/E4R Cr#kr+i#i) 

IF vr=O THEN p(tj)=FIl/2s30TO 250 
240 p(j)=ATN(i/r) 

250 IF p(j) + © THEN p{(j)=Fl+p(j) 
Q bn=oQ 

270 a=2#PFl#J/n 

#30 PRINT "j="s) 

270 FÜR i=0 TO n-1 

300 bn=FNfCi)#SIN (ati) +bn 

310 NEXT :i 

320 bn=£#bn/n 

330 y(j)=bn 

340 NEXT j 

350 MODE = 

360 LOCATE S:ZO:PRINT"f/fc=":k 

370 LOCATE 5:21 3: FRINT"m=":m 
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330 
590 
490 
410 
470 
450 
440 
450 
460 


100 


LOCATE S:22 : PFRINT'"'calcul sur “"snh 


FOR t=0 TO 400 

FOR J=O TO vh 

s=y(j)#t Cj)#SINIE#PIHj#t/600-p (j)) 
ss=ss+s 


NEXT 
FLOT 
ss=0 
NEXT 


j 
t-Z200+30+ss 


FNf(i)=SGN(SIN(Z*P1I#i/n)) 
nombre de points de calcul: 64 


Ô© 1.70141E+38 


FNfCIi)=SGN(SIN(Z#PI#i/n)) 
nombre de points de calcul: 44 


rapport f/fc 3: O.1 
O © © 


FNfCIi)=SGN(SIN(Z#PI#i/n)) 
nombre de points de calcul: 44 
rapport f/fc : O.1 

coeff. d amortissement: 0.2 


FNfCIi)=SGQN(SIN(Z#PI#i/n)) 
nombre de points de calcul: 44 
rapport f/fc : OQ.1 

coeff. d amortissement: 0.7 





f/fe= 6.1 
a= 8. 
calcul su 32 haraoniques 





3: "harmonique 


ENT C1) =SGNCESINCZ#FI# 1 /n)) 
nombre de points de calcul: 44 


rapport f/f5 3: O.1 


coeff. d amortissements 0.5 


f/fc= 0.1 
n= 6. 
calcul sur 32 harnoniques 


FNfCi)=SGNCSINCE#FI#i /n)) 
nombre de points de calcul: 44 
rapport f/fc % O.1 


coeff. d amortissement: 0.707 


JAY a 8.1 
a= 8. 
calcul sur 32 harnoniques 
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Transformée de Fourier rapide 


La transformation de Fourier est un outil très puissant dans l'étude 
des signaux de toute sorte. 


Le théorème de Fourier proprement dit peut faire correspondre à 
une fonction continue et périodique f (t) (de période T} le développe- 
ment : 


N= = N= © 
f{t) = Ao + 2  A,cosnwt + ÏZ  Bnsinnwt 
N=1 N=1 
où l’on a : 
1/7 
Ao = — f (t) dt 
T 
0 
= 
2 di 2 f (t) cos nwt dt 
B; = — f (t) sin not dt An = — 
TJ, T0 


Des expressions équivalentes sont les suivantes : 


n= 
e ft) = A + Z  aæ, cos (nwt + yh) 
n=1 
n=+ 
of(t) = Z  C,exp (nuit) 


n=- 


Cette dernière expression fait apparaître des fréquences «négatives », 
dont la présence se justifie par le fait que f (t) est réelle. Cela a pour 
conséquence l'existence de symétries, utiles pour la compréhension des 
phénomènes étudiés plus loin. 
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b 
e Ona a? = a2 +b2 n = — Arctg (=) 
n 
T/2 


e Ona C, -f f{t) .exp (— jnwt) dt = a, — jb, = ce 
T/2 


Ceci permet de faire correspondre au signal f (t) — périodique rap- 


pelons-le — une suite de couples de nombres : O0 CA 
1 ; 
2 œ 


qui sont les amplitudes des raies du spectre de f (t) ; ce spectre contient 
ainsi les pulsations O (le «continu »), 1 (le fondamental), 2, 3,..... | 


Que se passe-t-il si f (t) n’est pas périodique (soit si T + æ) ? On 
peut voir que les relations de définition des coefficients deviennent des 
intégrales et par conséquent des fonctions continues de w. 


On associe alors à f (t) la fonction F (w) définie par exemple par : 


+ 


F (w) { f{t).exp (— jwt) dt 


F (w) est la transformée de Fourier de f (t). 


L'échantillonnage d'un signal F (t) 


Considérons une fonction continue du temps f (t) et supposons que 
nous en notions les valeurs à des intervalles de temps égaux : Te. 


Au signal continu f (t), correspond maintenant un ensemble de 
valeurs numériques (f (0), f (Te), f (2Te) ...). 


On démontre que ces deux notions sont identiques si la plus grande 

Fe 

fréquence contenue dans le spectre de f (t) est inférieure à —— = — 
(Théorème de Shannon). 2Te 2 


Supposons maintenant f (tr) quelconque et tentons à l'aide d'un 
calculateur numérique d'évaluer la transformée de Fourier de ce signal 
f (t). C'est possible puisque nous possédons des nombres (les valeurs 
numériques de f (t) aux instants d'échantillonnage) et nous pouvons 
donc évaluer numériquement ces intégrales. 


Mais un problème surgit : il n’est pas possible d'évaluer ces intégra- 
les de moins l'infini à plus l'infini. Aussi est-t-on amené à limiter le 
nombre d'échantillons à un nombre fini N. 
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Cette limitation va introduire des erreurs : tout se passe comme si 
l'on multipliait la fonction f (t) par une «fenêtre», c'est-à-dire une 
fonction valant 1 de O à NT, = T et O0 partout ailleurs. Le spectre 
obtenu par une transformée de Fourier discrète (celle que nous tentons 
de réaliser) sera donc celui fourni par une «vraie» transformée, qui cor- 
respond au produit de f (t) par la fenêtre. On peut réaliser celle-ci pour 
savoir quel en sera le résultat. Celui-ci est simple : on reporte en chaque 
point de la «vraie» transformée de f (t) la transformée de la fenêtre et 
l'on effectue la somme des courbes obtenues (cette démarche s'appelle : 
convoluer les 2 fonctions). 


D'autre part, et ceci est toujours vrai, échantillonner un signal 
produit un spectre aussi par convolution : le spectre du signal original 
est reporté sur fe, 2 fe, ..., (on voit ici l'intérêt des fréquences « néga- 
tives»}) le spectre résultant étant la somme de tous ces spectres. Ainsi le 
spectre original n’est pas déformé (entre O et Fe/2) si effectivement la 
fonction étudiée ne contient pas de fréquences supérieures à Fe/2. 


Enfin, cherchant numériquement une fonction F (w) continue, il 
ne nous sera possible de n'obtenir que quelques valeurs de celles-ci. || 
est classique à l'heure actuelle d'effectuer une Fast Fourier Transform 
pour laquelle : 


— le nombre d'échantillons N est une puissance de 2 (4,8,16,32,64, ..). 


1 2 3 
— les fréquences pour lesquelles F (w) est calculée sont O, T' T' 7 





TT” en appelant T la durée de la fenêtre. 


N 
— la fréquence d'échantillonnage étant par conséquent T: le spectre 
obtenu sera toujours symétrique par rapport à FT : Seule la moitié 


des raies du spectre apporteront une information (dans la mesure où il 
n'y aura pas de recouvrement d’ailleurs, c’est-à-dire mélange des raies 
du spectre de f (t) avec les raies parasites créées par l'échantillonnage). 


Le problème de la fenêtre 


Tout se passerait fort bien si dans la fenêtre de durée T, il y avait un 
. : : : 1 2 3 
nombre entier de périodes des sinusoïdes correspondant à T' T' TT 
Malheureusement, ce ne peut toujours être le cas (signal non pério- 
dique, par exemple). 


Que se passe-t-il alors dans le cas contraire ? 


104 


Un retour à la méthode de convolution nous l’apprend : des raies 
parasites vont naître et fausser le diagnostic. Les mathématiciens propo- 
sent comme remèdes d'autres types de fenêtres dont les ondulations 
sont beaucoup plus faibles que pour la fenêtre rectangulaire, mais dont 
malheureusement la largeur du lobe est plus grande. 


Dans la fenêtre Dans le spectre 
1 
— 1 période de sinusoïde — 1 raie à T 
À : : a: 
— 2 périodes de sinusoïde — 1 raie à T 
; ; . ; 3 
— 3 périodes de sinusoïde — 1 raie à T 


On pourra constater sur les dessins l'influence de la fenêtre sur le 
spectre obtenu. 


Dans la figure 1, la fenêtre de largeur T contient un signal périodi- 
que, de fréquence fondamentale 1/T. Une seule raie apparaît alors dans 
le spectre, et c’est la bonne. 


Dans la figure 2, nous avons supposé que la fenêtre contenait 1,5 
périodes d'un signal. Ceci se traduirait normalement par une raie dans 
le spectre à 1.5/T. Mais la FFT ne peut mettre en évidence une telle raie ; 
elle ne calcule en effet que l'amplitude des raies présentes en 1/T,2/T, . 


L'utilisation de la fenêtre rectangulaire se traduit alors par l'appari- 
tion de fausses raies en 1/T, 2/T.. 


En utilisant un autre type de fenêtre, dont le spectre est plus « lourd » 
mais d'ondulations beaucoup plus faibles, il est possible, au détriment 
de la qualité du spectre d’un signal pur (figure 3), de réorganiser celui 
d'un signal de fréquence 1.5/T (figure 4). 


Cela est du au fait que la FFT ne sait pas que le signal a sa propre 
période ; elle considère que celle-ci est T dans tous les cas. 


Ainsi, le spectre observé n'est pas celui du signal proprement dit, 
mais d’un signal qui a la même forme que celui étudié, mais qui serait 
périodique de période T. 
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0 1/T 


1 
En continu : la raie à = (sinusoïde pure). 


En pointillé : le spectre de la fenêtre rec- 
tangulaire. 
La convolution de ces 2 fonctions ne 


2 3 
génère aucune raie à —,—,.… 
TT 





1,5 
Une raie à —— , non mise en évidence par 


la FFT génère des raies parasites lors de la 
convolution avec le spectre de la fenêtre. 





® Idem à (1) , mais fenêtre de Bilackmann. 


@ Idem à (2) , mais fenêtre de Blackmann. 





106 


nombre d’echantillons = 24 


20 FÜR i=1 TO mixr Ci) =SGNSINCS#PIHI HZ T/N)+SINCE#PIHIHIS, 
S/n)):xi(i)=O:NEXT 


ASFECT TEMPOREL 





aspect frequentiel 





nombre d’echantillons = ZE 
FENETRE LE BLACEMAN 


ASPECT TEMFÜUREL 
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aspect frequentiel 





49 FÜR i=0 TO nixr Ci) <SGNCSINCZEFIæi#1E,S/nm))exi CI EOENEXT 
nombre d'echantillons = 254 


ECT TEMPOREL 








aspect frequentiel 
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nombre d’echantillons = 254 


FENETRE DE ELACEMAN 


ASPECT TEMPOREL 





FENETRE LE ELACKMAN 


aspect frequentiel 
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60 FOR i=0 TO nixr(i)=SIN(2#PI#i#12.5/n) +SIN(Z#PI#i#2.5/n) 
®xi(i)=O:NEXT 


nombre d’echantillons = €4 


ASPECT TEMPOREL 


Il | | l | { | | | 


aspect frequentiel 
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romtre d'echantillons = 2546 


FENETRE LE EHLACEMAN 


ASFECT TEMPOREL 





aspect frequentiel 


nombre d'echantillons = 44 
FENETRE LE ELACEMAN 


ASPECT TEMPOREL 


gl 
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aspect frequentiel 


QD FÜR i=Q TO nixr Ci)=SIN(E#FIHi#12.5S/n)sxi Ci) OI NEXT 
nombre d'echantillons = 44 


AZFECT TEMFOREL 


We Le pi, 


aspect frequentiel 
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nombre d'echantillons = 44 
FENETRE DE EBLACKMAN 


ASFECT TEMPOREL 


L: an 
a m “nf | | DU ” 


aspect frequentiel 


60 FÜR i=0 TO nixrCi)=SIN(E*PI#i#12/n)txi Ci) =O:NEXT 


nombre d'echantillons = 64 


ASPECT TEMPFÜREL 


A dal " , 1 , j 1} | 
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aspect frequentiel 


nombre d’echantillons = 64 
FENETRE LE BLACKMAN 


ASPECT TEMPOREL 


_- “rl ph Li. 


aspect frequentiel 
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nombre d’echantillons = 64 


ASPECT TEMPOREL 


vo HU 


aspect frequentiel 


60 FOR i=0 TO nixr(i)=SGN(SIN(2#PI#i#2/n)):xi (i)æO:NEXT 


nombre d’echantillons = 16 


ASPECT TEMPOREL 
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SFECT TEMPOREL 


4 1 | 
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10 KREM it 
20 REM # # 
3 


QO REM *% fft * 
40 REM # # 
50 REM it 
40 RENM 
70 REM fft 


INPUT "nombre d’echantillons:2"n="in 


90 LIM xr(m)sDIM xi Cm) 


100 
110 
120 
130 
140 
150 
160 
170 
180 
190 
Z00 
210 
220 
230 
240 
250 


INFUT'"'puissance de Zi"s;nu 
ni=n/i 

FOR i=O TO mixr (i)=SIN(Z#PI#i4kZ/nm)sxi (i)=O:NEXT 
FRINT "fenetre de blackman ? (o/n)" 
INPUT cg 

IF cé="n" THEN 200 

FÜR izQO TO n 
q=2+#Pl#Ci-n/2) /n 
xr(i)=xr(i)#(0.42+0.S4C05S (q)+0.084C05 (Z#q) ) 
NEXT i 

MODE 2 

FÜR 1=0 TO N-1 

PLOT I#xSI1Z2/N: 200 +1O00%XR(I) 
DRKAW I#512/N: 200 

NEXT 

ni=nu-i 

k=Q 

FOR 1=i TO nu 

FÜR i=1i TO n2 
a=INT(kK/(2"n1)) 

b=nu 

GÜSUR 680 

p=br 

arg==#PI%xp/n 

c= COS (arg) 

Ss=$SIN(arg) 

tr= xrCk + nZ)kc+xi (K+n2)%+s 
tis=xi (ktn2)#c-xr (k+n2)*%s 
xr(k#n2)æxr (K)-tr 
xi(k+n2)=xi (K)-ti 
xr(k)=xr(k)+tr 

xiCKH)=xi (kK)+ti 

kK=k+1 

NEXT ii 

k=k+n2Z 

IF k<n-1 THEN GÜTO 280 

k=O 

ni=ni-i 

n2 =INT(n2/2) 

NEXT 1 

FOR k=O TO n-i 

a=ktbhænu 

GOSUB 680 


117 


530 
540 
550 
S60 
S70 
580 
590 
600 
610 
620 
630 
640 
650 
660 
670 
680 
690 
700 
710 
720 
730 
740 
750 
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i=br 

IF ik THEN tr=xr(k)itizxi (KD)3=xr(K)=xr(i) 
IF ik THEN xi(k)=xi(i)sxrCi)=trasxi(i)=ti 
NEXT K 

MOLE 2:CLS 

FOR i=0 TO n-1 
xr(i)=SGR(xr Ci) 72+xi (1) 72) 
IF xr(i) > MX THEN MX=xr(i) 
NEXT 

PLOT O:1:DRAW S12:1 

FOR x=0 TO n-1 

PLOT x#512/n:0 

DRAW x*512Z/n; 320%xr (x) /MX 
NEXT 

STOP 

ji=a:z=0 

FOR x=z1i TO b 

j2=INT(j1/2) 
z=z#2+(j1-j2+2) 

ji=j2 

NEXT x 

br=z 

RETURN 


Nombres aléatoires 


Beaucoup d'applications informatiques nécessitent l'emploi de nom- 
bres aléatoires, c'est-à-dire fournis par le hasard : jeux, simulations, … 
Tous les calculateurs, gros ou petits, disposent donc d'un générateur de 
tels nombres. 


Malheureusement, ces générateurs ne sont pas parfaits et ne peuvent 
l'être. Ils sont obligés d'utiliser un algorithme, c’est-à-dire une méthode 
mathématique, fournissant une suite de tels nombres. Ceux-ci ne sont 
donc pas réellement fournis par le hasard, mais obéissent à une loi : si 
celle-ci a été bien choisie, ils seront apparemment aléatoires (d’où 
l'appellation de « pseudo aléatoires »). 


L'Amstrad possède un tel générateur ; c'est la fonction RND 
(RANDOM), qui peut posséder un argument. 


— RND (X) avec X > O0 ou RND fournit une suite de nombres, aléatoi- 
res compris entre O (inclus) et 1 (exclus). Aprèsun RUN, la suite fournie 
débute toujours par le même nombre (0,140522114). C'est l’utilisation 
classique de cette fonction. 


— RND (0) fournit le dernier nombre généré par RND. 


— RND (X) avec X < 0 génère une nouvelle séquence aléatoire et fournit 
le premier nombre de cette série. 


Le programme 1 montre à l’aide d'un exemple le fonctionnement 
de ces instructions. 


En pratique, on peut avoir besoin d’une série de nombres aléatoires 
sur un intervalle différent de 0,1. (La simulation d’un jet de dés nécessite 
des nombres entiers compris entre 1 et 6). Pour cela, la transformation 
des nombres e fournis par RND se modifie aisément : soit a la borne 
inférieure de l'intervalle recherché et b sa borne supérieure. Alors : 


r =e(b —a) + a est un nombre aléatoire compris entre a (inclus) et 
b (exclus). 
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= [INT (RND * (b —a + 1) ) + a est un nombre aléatoire entier (si 
a l'est) compris entre a et b (bornes incluses). 


Ainsi, pour un jeu de dés : (a = 1,b = 6) :r = INT (RND * 6) + 1. 


L'exemple 1 montre ainsi la génération d'une suite réalisée avec 
RND, fournissant des nombres compris entre 1 et 10. La première 
colonne indique le numéro du tirage, la seconde le nombre correspon- 
dant et la troisième la moyenne des tirages. Celle-ci devrait valoir 5,5 
1+10 


soit 
2 


Programme 1 


10 FOR N=i TO 10 
£O PRINT RNO(1) 
30 NEXT 

35 FRINT 

40 FÜR N=i TÜ 10 
SO PRINT RNL(G) 
&Q NEXT 

6% PRINT 

70 FÜR N=1i TO 10 
&Q PRINT RNL(-1) 
70 NEXT 

9 PRINT 

100 FOR Nz 1 TO 10 
119 PRINT KRNL(1) 












1:0 NEXT 
0.271974005S5 
0.140 Fe de O.:7197404 
_ &7% QE O0.271740645 
Age Q0.2717404 


0.27174045 
Q.2717404 
0.27174065 
0.77174094 
0,.271740455 
0,271974042385 


O.S2SA1ZSEE 
O.021330127 
0.175135014 








0.630562254 
QO.43054Z2S4 
0.650562:54 
© a 








pe 491035 





3 
PRET 14, 





OO ASOSALTE ï 
O.43056L254 0.78951017%3 
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borne inferieure : 1 
borne superieure =: 10 

10 9.594671556 5.51575454 
£9 5.401230742 5.73377248 
30 4.8%9457:744 4,.0177656S 
49 1.42745654 5.8757S22 
50 S.1808817S S5.6442%027 
40 7.43551282 5.97205253 
79 7.44405S28 5.81507795 
Este) 2.07373057 5.705229729 
90 S.4957Z42E6 5.771301 
100 7.10377051 5.84357446 
119 7.30424957 2.731:4674 
120 4.10745742 5.74784555 
130 Z.5730447 S5.73884501 
140 9.79324457 5.77024574 
150 7.14415371 5.81574419 
160 8.472582545 5.31750322 
170 Z.38825701 5.23954354 
130 1.71401477 5.313%7739 
190 3.01:03:09 5.77044382 
Z90 2.3411512 5,.77347373 
10 3.120975435 5.723223777 
Z20 2.16051473 5.49254442£ 
Z30 4.0259L4544 5.65741453 
Z40 3.54755336 5.47114233 
220 4.41330752 5.48473544 
Z40 5.47259977 S.47561453 
279 3.07744425 S.6444553748 
ZS0 2.7971974473 5.79177174 
270 1.401145 5.47444178& 
300 4.27308257 5.647465658 
310 5.7175112 2.7321342 
320 7.27182538 5.7600079S 
330 Z.45056251 5.744457 
340 5.7:0457274 5.74593401 
330 D.42197:474 3.750577S6 
340 7.333224L04 S.73350317 
370 2.74177747 5.74710013 
330 1.23347974 5.423274149 


Le programme 2 donne un exemple de générateur de nombres aléa- 
toires, s possédant les mêmes propriétés que RND. L'algorithme figure 
en lignes 100 à 120. Il nécessite une sources (initialisation de la série). 
L'exemple 2 montre l’utilisation de ce programme. 
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Programme 2 


10 REM ##ketetetoteote sheet cheminee ee RAGE 
O REM * # 
30 REM * DISTRIEUTION UNIFORME SUR AB *# 
A4Q REM * # 
DO REM #steteitestonitestemeneieteoetesetenpetetese teste Meet ee 
&O REM 

70 INPUT "borne inferieure : 

30 INPUT "borne superieure :":b 
9Q INPUT “ source :";s 

100 e=(s+24:98+99991)/1997017 
119 e=e-INT(e) 

120 s=e+197017 

130 r=e*(b-a)+a 

140 PRINT vr ; 

145 m=mtrin=n+i 

146 PRINT m/n 

150 GÜTI 100 

aleal 1 
























torne inferieure : 1 

borne superieure : 10 

source s 1 

19 5.5223546444 

£0 4.7407074 

30 3.7333437 

40 4.53052521 

50 2.642:149371 

(21e) 7.295110786 

79 7.1365E2574 

80 3.35451354 

70 b.LESSZISS1 

100 7.Z4133301 5 

119 $.7%747078 : & 

120 3.53147146 5.756 + 
130 3.51079376 5.787 4 
140 6.11774%707 5.841517537 
150 3.434536c74 5.753173745 
140 2.446433454 2.672335315 
170 4.4545232S LE 9034 
150 Z.:335:432 SLEY2z 
170 2.277567 ‘43 & 
Z00 3.072544 3 
£10 4.147348 £ 
Z20 2.517754 

Z230 1.07093043 

#40 7.5717571S 

E5O 7.74557114 

Z40 5.74035149 

270 £.51104449 
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tants. On dit alors que l’on a affaire à une distribution uniforme. 
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5.4404528 
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On peut visualiser la qualité de tels générateurs à l’aide du program- 
me 3. L'intervalle (0,1) est découpé en M intervalles égaux ét après 
chaque tirage, le nombre fourni est placé dans l'intervalle auquel il 
appartient. Au bout d'un certain temps, l'aspect de l'enveloppe de la 
courbe devrait être à peu près horizontal, sans trous ni bosses impor- 


Programme 3 

EM either tee ME cr he 4e ee cd 4e fe 24 de 2e ce 24e 24 ee ee ge 
EM * # 
EM # ETULE DE LA DISTRIBUTION UNIFOÜRME * 
EM *% * 
EM tte edenhecte ste este che ohe see se pete 28e ee ee ee 
EM 

CLS 


REM distribution uniforme 
INPUT "nombre d’'intervalles : ":M 
INPUT "nombre d’iterations : 
LIM am) 
FOR i= 1 TO IC 
sn=RNE (1) 
su=su+sn 
av=i 

mn=mn+sn/av 

ss=ss+sn+sn 
va=vatsn#sn/(av#av) 
sn=(sn/av)*m 

FÜR j=0 TO m-1 

IF INT(sn)= j THEN 150 

NEXT j 

GSOTO 40 

a(j)=atj)+i 
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170 
200 





2350 


40 
nl €) 
z40 
270 


230 
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IF atj) > 150 THEN 170 
icr=INT(6é40/m) 

FLOT j#icrs:a(j):DURAUR icr:0 
NEXT :i 

i=i-i 
mn=mn/iisd=SQR((va-mn#mn/i)/i) 
FRINT'MOYENNE =": 

FRINT USING "#.###"; mn 

FRINT "LEVIATIUON STANDARD ="3 
FRINT LISING"#.###":sd 


nombre d’intervalles : 7 19 

nonbre d’iterations :7 10060 
1=100 

nombre d’intervalles : 7 10 

nombre d’iterations :7 1000 
1= 300 





nombre d’intervalles 
nonbre d’iteration 


3 7 10 
:2 10660 
1= 900 


Cependant, il n'existe pas que ce type de distribution dans la nature. 
Si la distribution normale convient pour les jeux de hasard (dés, loto, ….), 
elle ne s'applique pas à d'autres classes de phénomènes. Ainsi, la distri- 
bution des tailles d'individus dans une population donnée, si elle dépend 
du hasard également, n'est pas pour autant uniforme : si en moyenne 
elle vaut 1,75 m, il y a peu de chances de trouver des individus de 2,50 m 
en 0,50 m.. 


De telles distributions sont dites normales et sont caractérisées par 
la valeur moyenne mn de la variable x; et par l'écart type (en déviation 
standard) sd. 


1 
mn =—{(x; +xit...+x,) 
n 


sd = SOR ((n, — mn}? + (x2 — mn)? +...+(x, — mn?) 


sd renseigne sur la dispersion de la variable autour de la valeur 
moyenne mn. 


Il est remarquable que l'on puisse générer de telles suites à partir 
d'une distribution uniforme. Nous le verrons plus loin. 


Quel peut être l'intérêt pour un programmeur d'avoir à sadisposition 
une distribution normale ? I! peut être le suivant. Considérons un mobile 
qui doit se déplacer sur l'écran du moniteur, dans une direction déter- 
minée. En utilisant des petits déplacements obéissant à la loi normale, 
ce mobile peut donner l'impression de se déplacer de façon aléatoire, 
alors qu'en fait il le fera — sur une longue distance — dans la direction 
choisie. 


Revenons à ce type de distribution et tentons d'en réaliser une 
grâce à la fonction RND. Pour cela, prenons av nombres aléatoires suc- 
cessifs dont nous faisons la moyenne mn, puis de nouveau les av nom- 
bres aléatoires suivants, dont nous faisons également la moyenne mn 
et ainsi de suite. Nous disposons ainsi d'une nouvelle suite de nombres 
qui n'obéissent plus à une distribution uniforme. En effet, chacun 
d'entre eux n'est plus capable de prendre n'importe quelle valeur com- 
prise entre O et 1, mais plutôt une valeur proche de 0,5, valeur vers 
laquelle tend mn si av est très grand. 


Nous créons ainsi une distribution dite «normale», dont la valeur 
moyenne est 0,5 et dont l'écart type est réglable par av. Le program- 
me 4 permet d'examiner dans le temps l’évolution d'une population 
régie par une telle loi, (voir l'exemple 4). 
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R 
R 


R 
R 
R 
R 


[CRE 


Ho O0 à 


Programme 4 


EM ete ske she tee fe ee CEE EE NA A SRE MEN 
EM * # 

EM * DISTRIBUTION NORMALE * 

EM * + 

EM texte tete she be he 29e Re RENE RARE SES A SENS SEE MERE 

EM 

CLS 

INPUT'"nombre d’intervalles : ":M 
INPUT "nombre d'iterations :"3:1C 


INPUT " moyenne effectuee sur : ?";av 
LIM am) 
FOR i= 1 TO IC 
sn=0 
FÜR t= 1 TO av 
sn=sn+RND(1) 
NEXT t 
su=su+sn 
mn=mn+sn/av 
ss=ss+sn*xen 
va=va+tsn+sn/(av+av) 


160 sn=(sn/av)*m 


FOR j=0 TO m-i 
IF INT(sn)= j THEN 210 
NEXT j 


ZOO GOTO 270 


a(j)=a(j)+i 


NI 


ZE0 1F a(j) > 150 THEN Z00 


icr=INT(640/m) 


240 FLOT j+icrra(j):DURAWR icr:0 
250 NEXT i 

260 i=i-i 

270 mn=mn/i:sd=SGR((va-mn#mn#i)/ 1 
230 PFRINT'MOYENNE ="; 

90 PRINT USING "#.###"3: mn 
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QG ON 
= © 
© © 


PRINT "LDEVIATION STANDARD ="; 
PRINT USING'"#.###":sd 


MOYENNE= 0.49795104 
D'EVIATION= 0.1736019E3 
nombre d’intervalles : 10 
nombre d’iterations : 1000 


moyenne effectuee sur : 3 


127 


Le programme 5 permet d'étudier la répartition des nombres aléa- 
toires en fonction de av. Comme prévu, l'écart type diminue quand av 
augmente. 


Il existe d'autres algorithmes permettant d'obtenir une distribution 
normale à partir d'une distribution uniforme. Par exemple, si n; et 
n;,, sont 2 nombres aléatoires successifs (générés par RND), alors la 
suite : 


sn = V —2.Log{(n;) * cos (2.7.n;,:) 


obéit à une loi normale, de moyenne 0 et d'écart type 1 : la distribution 
est dite centrée et réduite. De façon plus générale, on peut centrer une 
distribution normale en soustrayant la valeur moyenne à tous ses élé- 
ments et on peut la réduire en divisant toutes les valeurs de la variable 
par l'écart type. 


Inversement, partant d'une distribution centrée réduite (x;), on 
peut obtenir une distribution normale (x;) de moyenne M et d'écart 
type ET par la relation : 


x; = xj.ET+M 


Enfin, pour ceux que cela intéresse, rappelons que la fonction de 
densité f (x) d'une distribution normale est donnée par la loi : 


1 
f(x) = — . exp (— (x; —M)?/(2ET?)) 


V27 ET 


Programme 5 


10 REM  :##esteite st etes se se ee tee 


Z0 REM * “ 
SO REM * loi normale * 
40 REM * # 
SO REM teinte tk 
6Q REM 


79 INK O,O:INKH 1:26 
50 PAPER Q : FEN 1 


90 MODE Z 

100 CLS 

110 REM distribution normale 

120 INPUT'"nombre d’intervalles : ":M 


130 INPUT "nombre d’iterations :":1C 
140 INPUT "moyenne effectuee sur :":av 
150 DIM am) 

160 FÜR i= 1 TG IC 

170 sn=0 

180 FÜR t= 1 TÜ av 
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500 
510 
SzQ 
580 


540 


sn=sr+RNEL (1) 

NEXT t+t 

su=su+sn 

mi=mnh+sn/av 

Ss=ss+sn+#sn 

va=va+sn#sn/ (av#av) 
sn=(sn/av)#+m 

FÜR j=0 TO m-i 

IF INT(sn)= j THEN 300 

NEXT j 

GOT 340 

a(j)=a(j)+i 

IF atj) * 150 THEN 270 
icr=INT(é4Q/m) 

FLOT j#icr:a(j):DRAWR icr:0 
NEXT :;i 

i=i-i 
mi=mn/i:sd=SQR((va-mn#mnxi)/i) 


FRINT'MOYENNE =": 

PRINT USING "#.###"5 mn 

PRINT "DEVIATION STANLDARL =": 
FRINT USING"#.###":sd 
a=m*su/ay 

b=ss+m+m/av/av 
ss=SOR((b-A#A/i)/i) 

si=sd 


d=ss+k£aR (2*PI) 

sz=si+si 

RKEM i/d est le max. de la courbe 
PLOT O,i/d:DRAWR 600: 0 

FOR x= OO TO 1 STEF 1/640 

numsEXF (—(x-mn)7"2z/ (Z#s2)) 
f=x4+6@40 

g=num#i/d 

PLOT f,g-:1 

NEXT 
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130 


nombre d'intervalles : 19 


nombre d'iterations * 1000 


moyenne effectuee œur ? 7 


MOYENNE -6,498 
DEVIATION STANDARD -8.198 





ombre d'inutervalles : 10 
nombre d'iterations ? 1000 


moyenne efferituee sur 3 à 


MOYENNE 6.500 
DEUIATION STANDARD -6.166 





nombre d'intervalles * 10 


nombre d'iterations : 1900 


moyenne effectuee Sur : 4 





nombre d'intervalles : 10 
nombre d'iterations : 1000 


moyenne effectuee sur ?: © 


MOYENNE =9.584 
DEVIATION STANDARD -8.134 
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nombre d'intervalles : 19 
nombre d'iterations : 1000 


moyenne effesctuee sur =: € 


MOYENNE =8.493 
DEUIATION STANDARD -8.114 





nombre d'intervalles : 19 
nombre d’'iterations : 1000 


moyenne effectuee sur : 7 


MOYENNE =8.581 
DEUIATION STANDARD -8.118 





nombre d’'intervalles : 10 


mombre d’'iterations : 1900 
moyenne effectuee sur ? & 
MOYENNE =-0.497 UN 
DEVIATION STANDARD -8.106 / \ 





Terminons par une petite application de toutes ces notions : le calcul 


de x. Considérons un carré de côté a ; sa surface vaut a? . De même, soit 
2 


a a 
un cercle de rayon 3 ; Sa surface vaut 7 FE Conclusion : le rapport 
us ; 
des surfaces vaut Fi Plaçons un point au hasard dans le carré ; la proba- 
: ne T 
bilité qu'il tombe dans le cercle est évidemment PE On peut donc éva- 


luer 7 par le rapport : 


” nombre de points dans le cercle 
T = ————————————————————— 
nombre de points dans le carré 


D'où le programme 6 et l'exemple 6. 
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10 
Z0 


70 


100 
110 
120 
130 
140 
150 
160 
170 
1aQ 
190 
200 
210 
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Programme 6 


REM _#24ete tete etes tee te 4e ee 


REM * *# 
REM + calcul de pi * 
REM * * 
REM _ #24:4658 4e 4e AO HO MONTE 
REM 

MOLE 2 

PLOT 170; SO:URAUWR O:,30Q 


DRAWR 300,0: DRAWR O:-300:DRAWR -300; © 
FOR n= © TO Z#PFI STEF P1/300 

PLOT 320+150*#C0S (n) ; ZOO+150#SIN (mn) 
NEXT n 

FÜR i=0 TO 100000 

x=RNEO (1) 

ÿY=RND (1) 

IF x%x +y#y < 1 THEN p=p+i 

PLOT 170+x#300; 50+y*#300 

IF i/100 <> INT(i/100) THEN GOTA Z0Q0 
LOCATE 1,1: PRINT"'pi=":p#4/i 

LOCATE 1» 2:PRINT"i=";3 ; 

NEXT i 





i= 700 


i= 1000 


pi= 3.08 07908 
i= 1980 


i= 1100 
i= 1200 
i= 1300 
i= 1400 
iz 1500 
iz 1600 
i= 1700 
iz 1€00 
i= 1700 
i= 2000 


pi= 3.1 954978 
i= 2806 


3.072664LE7 


= 3.03 


= 4.04 


4.05333333 


3.02857143 


3.015 


3.07672305 


3.08571429 


5.09706666&7 


EN] 
ù 

[e] 
cc 


47 OS 


= 3.095556 


3.07684%11 


3.1 
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He: /2109 pi= 3.11428S71 


i= 2#00 pi= 3.11434364 
i= 2500 pi= 3.11652174 
i= Z400 pi= 3.12333333 
i= Z200 pi= 3.156 
i= Ze00 pi= 3.14925077 
i= £700 pi= 3.15407407 
i= Z&00 pis 3.14 
i= 2900 pi= 3.13655172 
i= 3000 pi= 3.134 


pi= 3.136 119 
i= 3008 
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Valeur moyenne 
Valeur quadratique moyenne 


Rappelons tout d'abord les définitions de la valeur moyenne et de 
la valeur quadratique moyenne d'une fonction y (x) sur l'intervalle 


Xar Xb : 
a b *b 


1 
valeur moyenne : MO = ———, y (x). dx 
(xb = Xa) x 


Xb 
valeur quadratique moyenne : QM = sn . {y (x))2. dx 
(Xp — Xa) F 


ce qui implique : OM = valeur moyenne du carré de y. 


Dans ce programme l'ordinateur calcule la valeur moyenne et la 
quadratique moyenne de la fonction y (x) proposée par l'utilisateur 
sur un intervalle de son choix. En outre il représente graphiquement 
y (x) et (y (x))? et met en évidence leur valeur moyenne. 


Vous pouvez éviter de recopier les lignes 10 à 140 qui ne servent 
qu'à la présentation. 


Déroulement du programme 


L'utilisateur donne la fonction y (x) en ligne 1000. Puis il lance le 
programme en faisant RUN 200. Il indique alors l'intervalle de valeur de 
x étudié. 


L'ordinateur calcule y (x) pour 600 valeurs de x dans cet intervalle. 
Le nombre 600 correspond au nombre de pixels sur la largeur de l'écran 
(compte tenu des marges). |! calcule une valeur approchée de l'intégrale 
en faisant la somme des surfaces de 600 rectangles. 1! détermine alors la 
valeur moyenne. || procède enfin de même avec la fonction (y (x))?. 
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Tout en calculant les valeurs de y, il compare chaque valeur à la 
précédente afin de connaître les valeurs extrêmales de y sur l'intervalle 
considéré. ll en déduit une échelle de représentation en ordonnée pour 
tracer la courbe sur l'écran. Sont également tracées les droites indiquant 
les valeurs moyennes de y et de y?. 


Après avoir observé y (x) et sa valeur moyenne, on tapera sur la 
flèche pour obtenir la valeur quadratique moyenne. 


Détail du programme 


Lignes) 
10 préambule 
210 intervalle de x considéré 
300 boucle de calcul 
1000 expression de y (x) 
320 mise à jour du maximum et du minimum 


340 mise à jour de l'intégrale de y 

360 mise à jour de l'intégrale de y 

390 calcul des valeurs moyenne et quadratique moyenne 

410 détermination du domaine de valeurs de y 

510 on étend ce domaine à zéro pour représenter l'axe des 
abscisses sur l'écran 

540 graphe de y (x) 

640 graphe de (y (x))? 


Liste des variables 


x variable 
Xa,Xb bornes du domaine de x 
Xe largeur de ce domaine 


ex échelle en x 
dx accroissement de x entre deux calculs 
p,m valeurs maximale et minimale de y 

X 


i - f y (x). dx 


x 


-f OA de 


a 
z =" ÿ? 
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mo 
qm 
pp 
ey 


XX 


valeur moyenne de y 

valeur quadratique moyenne de y 

valeur maximale de y? 

échelle en ordonnée pour la représentation de y 
décalage en ordonnée pour placer l'axe des abscisses 
abscisse d’un point sur l'écran 


lflustrations 


de valeur de x :x, = 0, x, = 2.PI, soit 6,283185. 


Comparons les résultats obtenus avec les résultats théoriques : 


Figure 1 : 


y = sin (x) [[ MO = 0 par raison de symétrie 
0 = 1//2 
valeur mo enne = 9.74825E-09 









fi HN 


Fig. 1a 





quadratique moysehne = 6. 70716679 


Toutes les illustrations ont été réalisées en prenant pour intervalle 
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Figure 2 : 
= {sin (x)|[[ MO = 2/PI 
QM = 1/2 (on retrouve la valeur précédente car 
le signe disparaît dans l'élévation de y au carré). 


| A li 


valeur (i = -63661397 

















fn 
y 


Fig. 2 


























4 | 
































Figure 3 : 
y = sin (x) si x est inférieur à PI sinon y = 0 
= 1/PI (il n’y a qu'une seule arche, la surface 
sous la courbe est donc deux fois plus petite). 
OM = 1/2 (la surface sous la courbe y est deux fois 
plus petite donc le carré de la valeur quadratique 
moyenne est deux fois plus petit). 


valeur mouenne = 60.31830699 





| LL 


Fig. 3a 
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qauadratique moyenne = 60.500000011 


i À 


F \ ü AL 


Figure 4 : 
y = 2 + sin (x). 

MO = 2 (car la valeur moyenne du sinus est nulle). 
y? = 4+4,sin (x) + (sin (x)}? 

e la valeur moyenne de 4 est 4 

e la valeur moyenne de 4.sin (x) est nulle 

e la valeur quadratique moyenne du sinus est égale à 
1//2 donc la valeur moyenne de son carré est : 
(1/2)? = 1/2. 

e finalement OM = ;/4,5. 


valeur moyenne = 


ail (l D. . | 
+ 4 |. 

































































Fig. da 
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2.12132036 


1 





Programme 


10 CLS 

20 LOCATE 12,6: PRINT'"valeur moyenne" 
SO LOCATE 8,8:PRINT'"Valeur quadratique moyenne" 
40 LOCATE 1,22:PRINT'"pour faire avancer le program 
me, taper "CHR$(243) 

50 LOCATE 39,25:PRINT CHR$ (245) 

60 q=INKEY(1):1F q<>0 THEN 60 

70 CLS 

100 PRINT:PRINT'"donner la fonction y(x)" 
110 PRINT:PRINT'"en ligne 1000" 

120 PRINT:PRINT'"puis faire RUN 200" 
130 PRINT:PRINT:LIST 1000 

140 END 

200 CLS 

210 PRINT:PRINT'"valeur de x initiale " 
220 INPUT xa 

230 PRINT:PRINT'"valeur de x finale " 
240 INPUT xb 

2350 xc=2xb-xa 

260 ex=600/xc:dx=1/ex 

270 PRINT:PRINT'"calcul en cours" 

280 REM 

290 i=0:j=0:x=xa: GOSUB 1000:p=y:m=y 
300 FOR x=xa TO xb STEP dx 

310 GOSUB 1000 

320 IF y>p THEN p=y 

330 IF y<m THEN m=y 

340 i=i+dx*y 

330 z=y+y 

S60 j=j+dx#z 

370 NEXT x 

380 REM 

390 mo=i/xc 

400 qm=SGR(j/xc) 
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410 
420 
500 
510 
520 
330 
540 
990 
560 
570 
580 
590 
600 
610 
620 
650 
640 
659 
660 
670 
680 
690 
700 
710 
720 


pp= (MAX (ABS (p) ,ABS (m) ) ) 2 

REM 

CLS 

IF p<O THEN p=0 

IF m>0 THEN m=0 
ey=560/(p-m):k=ey*m 

MOVE 20,20-k:DRAWR 600,0 

FOR x=xa TO xb STEP dx#4 

GOSUB 1000 

XM=20+(x-xa) ex 

MOVE xx,20+motey-k:DRAWR O0, (ÿy-mo) #eyÿ 
NEXT x 

PRINT:PRINT'"valeur moyenne ="mo 
LOCATE 39,25:PRINT CHR$# (245) 
q=INKEY(1):1F q<>0 THEN 620 

CLS 

ey=560/pp 

MOVE 20,20:DRAWR 600,0 

FOR x=xa TO xb STEP dx*x4 

GOSUB 1000 

MOVE 20+(x-xa) #ex,20+qm*qmtey 
DRAWR O0, (y#y-qm*qm) *ey 

NEXT x 

PRINT: FRINT'"quadratique moyenne ="qm 
END 


1100 RETURN 
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Représentation d'une surface 
Z (x, y) 


Après avoir appris à dessiner des courbes, voyons ce qu'il en est des 
surfaces. Pour bien en concevoir la forme nous allons les représenter en 
perspective. On peut les suggérer par une distribution de points, mais 
l'interprétation de cette image est parfois difficile. Aussi nous allons 
tracer des courbes sur la surface. On peut tracer des courbes à x cons- 
tant, à y constant, ou à z constant. Cette dernière idée est difficile à 
mettre en œuvre. ici nous avons choisi de tracer des courbes à x constant. 


Si nous marquons tous les points ou toutes les courbes, l'image est 
confuse. Aussi nous allons considérer que la surface est opaque et ne 
marquer que ce qui est Vu ; ceci est un problème délicat. 


Premier essai 


Dans un premier temps, nous regardons la surface par le dessus 
(0 < ah < 90°) et la fonction z est telle que pour un couple x, y il n'ya 
qu'une seule valeur de z. Vous recopiez pour cela le programme jusqu'en 
ligne 300 en supprimant toutes les lignes dont le numéro ne se termine 
pas par zéro. 


Le cadrage a été choisi de façon suivante : x et y varient entre 
— 200 et + 200, z varie entre O et 200. Mais il n’est pas interdit d’envi- 
sager pour z des valeurs négatives. 


La direction du regard est définie par les angles 4h et av donnés en 
degrés suivant la figure 1. 
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Fig. 1. 


Envisageons d'abord la rotation d'angle ah. Le fait de ne pas être 
face au trièdre Ox, y, z est gênant, repérons-nous donc dans le trièdre 
O x,,ÿ:, 21. Le point de coordonnées x, y, z dans le premier trièdre a 


maintenant pour coordonnées x, y1,2z, avec (figure 2 ) : 


x, = x. cos (ah) — y . sin (ah) 
Y1 = x. sin (ah) + y .cos (ah) 


Z, —=2Z 





Fig. 2. 
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Tournons maintenant de av et considérons le nouveau trièdr 
Ox;,,ÿ2,2 avec : 


X3=X,.cos (av) +z, . sin (av) 
Yi = Yi 
Z2 = Z1.cos (av) — x, . sin (av) 


Le travail est fait pour chaque point de la surface en ligne 1300 à 
1330. 


Enfin la surface est vue à travers la fenêtre de l'écran. Les coordon- 
nées du point sur l'écran xe’ et ye sont calculées en litne 190 et 200. 
L'échelle e et les décalages 320 et 200 sont à choisir suivant la surface 
étudiée et l'angle sous lequel on la regarde. 


Ne marquer que ce qui est vu 


Le volume compris entre z = 0 et z (x, y) est découpé en tranches. 
Considérons la première tranche dont l’arête dessine sur l'écran (dans 
les coordonnées xe, ye) une courbe zz (ye). Dessinons maintenant le 
faîte de la seconde tranche. A l'abscisse ye, le point sur le faîte est à 
l’ordonnée ze, si ze est supérieur à zz (ye) le point est vu donc il est 
marqué, sinon il ne l'est pas. Si le point est marqué, ze est la nouvelle 
valeur de zz (ye) (figure 3). Ceci est réalisé en lignes 210 et 230. Les 
points sont reliés entre eux en ligne 240. Un point étant relié au précé- 
dent s'ils sont tous les deux visibles. 


Y [— courbe zz (ye) 






2e tranche 


1ère tranche 


Ye 


Fig. 3. 
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Nous pouvons calculer z pour toutes les valeurs entières de x et de 
y soit 400 X 400 = 60000 points. La haute définition de l'image 
d'Amstrad serait alors pleinement exploitée mais le programme serait 
très lent. Aussi nous avons choisi de prendre les valeurs de x de 15 en 15 
et les valeurs de ye de 4 en 4 (lignes 190 et 220). Aussi nous avons pris 
les valeurs de y de 4 en 4 en ligne 160. La haute définition n'est pas 
inutilisée pour autant puisque les points sont joints par des «droites » 
de haute définition. 


En ligne 40 nous vous proposons un exemple. Avec e = 1,ah = 30, 
av = 30, on obtient la figure 4. 


rt 


Br ; 
1j 





Fig. 4. 


Voir le dessous 


Si vous essayez des valeurs de av négatives dans le programme précé- 
dent, vous constaterez que vous ne pouvez pas voir le dessous de la 
surface, comme si elle était posée sur un socle opaque. Il faut donc faire 
un programme parallèle pour traiter le dessous. 


Soit t (x, y) l'équation du dessous, vous doublez le programme pré- 
cédent jusqu'en ligne 300 en recopiant les lignes manquantes. C'est 
lorsque te est supérieur à tt que le point est caché (figure 5). 





Fig. 5. 
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Dans un premier temps, regardons le dessous de la surface précé- 
dente : en ligne 45 vous placez pour cela : t = z. Avece = 1,ah = 30, 
av = — 15, vous obtenez la figure 6. 


2 


| l 
ie ù hi "il . 















3 n i Li il} su 0 il 
!, ou à u ! di « k Û n 
Qt “ sl 





Fig. 6. 


Mais vous pouvez aussi représenter un volume en choisissant une 
surface t différente de la surface z. 


Si vous placez en ligne 45 : 
t = 30 — 10 000/(70 + 1,3.MAX(ABS{(x),ABS{y))) 


avec les mêmes angles que précédemment, vous obtenez la figure 7. 





Les niveaux 


Reprenons la fonction z (x, y} et représentons-là vue de dessus à la 
manière d’une carte de géographie : plus z est grand, plus la teinte est 
claire. Pour cela vous recopiez la fin du programme et vous lancez le 
programme par un RUN 330 si vous ne désirez voir que cette dernière 
partie. 


Le terrain x, y de 400.400 est découpé en rectangles de 8 sur 16 
pixels c'est-à-dire la dimension d’un caractère en mode 2. La valeur de 
z sur ce petit rectangle est comparée à une échelle de 17 valeurs corres- 
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pondant à 17 caractères plus ou moins pleins. On aurait pu définir de 
nouveaux caractères avec une répartition de points au hasard dans le 
petit rectangle et une échelle de 8.16 = 128 valeurs plutôt que d'utiliser 
les caractères «tout faits» de l’Amstrad. 


On retrouve ainsi en figure 8 les trois bosses de la surface z (x, y). 


calcul en cours 





RZZALS : “bts . 
azit NÉS 
NATUITÉR ATUTES 
.“={fffiz: tzifff (= 
./ilfffl=: :=[fft0i/ 
NAT ETUUR 
./=[fflis:. { [ +izilff{=* 
SALES CORNE MDEE 
UN MEET ) } RAA à 
Fig. 8 
Organisation générale du programme 
Ligne(s) 
40 à 50 définition des surfaces 
60 à 130 introduction des angles 
140 à 300 dessin des surfaces 
170 calcul dezett 
180 envoi à 1300 matrice de rotation 
190,200 coordonnées du point z sur l'écran 
210 est-il vu ? 
220 il est marqué 
240 il est joint au précédent 
265 et 262 à 266 même travail pour le point t 
330 représentation par niveaux 
370 calcul de z 
390 liste des caractères 
400 adaptation de z à l'échelle des caractères 


410 à 440 recherche du caractère 
450 affichage 


Liste des variables 


e échelle 
ah,av angles de vue 
sh, ch,sv,cv lignes trigonométriques de ces angles 
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X,y,2,t coordonnées des points « dans la réalité » 

si d = 1 le point de dessus est vu 

sic = 1 le point de dessous est vu 

X1,V1,Z1,X2,V2,22 intermédiaires de calcul 

ye, ze,te coordonnées des points sur l'écran 

zz (ye) valeur maximale précédente de ze pour l'abscisse ye 
tt (ye) valeur minimale précédente de te pour l’abscisse ye 
zz,tt Valeurs précédentes de ze et te 

w correspond à z dans l'échelle des niveaux 

ww, n intermédiaires de calcul 


Programme 


10 REM'"surface z(x,y)" 

20 REM 

30 CLS:GOTO 60 

40 z=200% ((C0S (y#0. 024) ) 72) # (COS (y#0. 008) ) # ( (COS (x 
#0.008) ) 72) # (COS (x#0. 008) ) 

45 t=z 

SO RETURN 

60 DEG 

70 PRINT: INPUT'"'e=":e 

80 PRINT:INPUT"ah (-90< <90 )="3;ah 

90 PRINT: INPUT"av ( O< <90 )=";av 

100 Sh=SIN(ah):ch=C0S (ah) 

110 Ssv=SIN(av):cv=COS (av) 

120 RAD:CLS 

130 DIM 2z (160) 

135 DIM tt (160) 

136 FOR n=0 TO 160: tt (n)=400 

137 NEXT n 

140 FOR x=200 TO -200 STEP -15 

150 d=0 

153 c=0 

160 FOR y=-200 TO 200 STEP 4 

170 GOSUB 40 

180 GOSUR 1300 

190 ye=INT((320+e#y2)/4) 
200 ze=200+e*z2 

295 te=200+ext2 

210 IF ze<zz(ÿe) THEN GOTO 
220 PLOT 4#*ÿe,ze 

230 zz(ye)=ze 

240 IF d=1i THEN DRAWR ye-yy,ze-zz 
230 d=i 

260 zz=1e:yy=ÿye 

262 IF te>tt(ye) THEN GOTO 275 
264 PLOT 4#ye,te:tt (ye)=te 

266 IF c=1 THEN DRAWR ye-yy,te-tt 
268 c=1 

270 REM 


N 
œ 
el 


150 


275 
280 
290 
300 
310 
320 
330 
340 
350 
360 
370 
380 
3970 


5233 


400 


430 
450 


480 
900 
910 


1300 
1320 
1325 
1350 


tt=te 

NEXT y 

NEXT x 

REM 

LOCATE 539,25:FRINT CHR$(243) 
qg=INKEY#(1):1F q<>0 THEN GOTO 320 
MODE 2:PRINT'"calcul en cours" 
TAG 

FOR x=-180 TO 180 STEP 16 

FOR y=-200 TO 200 STEF 8 
GOSUB 40 

MOVE 330+y,200-x 

DATA 52,46,96,58,192,61,105,91,102,56,71,68,66 
1135,225,1435 
W=INT((z+13)/12) 

FOR n=1i TO w 

READ ww 

NEXT n 
PRINT CHRS$ (ww) ; 

RESTORE 

NEXT y 

NEXT x 
TAGOFF 

CALL &A000 
END 

REM 

xi=x#ch-y#sh 

tist 

yi=x#sh+y#ch 


1340 zi=z 
1350 x2=xi#cv+zl#sv 
1355 t2=tilxcv-xi#xsv 


1360 


y2=yi 


1370 z2=z1#cv-x1#sSVv 


1380 


RETURN 
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Représentation d'une surface 
en coordonnées cylindriques 


Nous nous proposons maintenant de représenter une surface définie 
par une équation en coordonnées cylindriques. Le volume enfermé par 
la surface est d'abord découpé en tranches (comme une orange), puis 
découpé en disques. 


Compromis 


La technique utilisée pour représenter la surface z (x, y) est diffi- 
cilement utilisable ici. Nous devons nous orienter vers une autre métho- 
de qui nécessite la mise en mémoire de l’état des points de l'écran. Ceci 
utilisant une mémoire très importante, nous sommes amenés à travailler 
en basse résolution ce qui procure d'autre part l'avantage d'un travail 
rapide au détriment de la qualité de l’image. C'est en tout cas un bon 
compromis pour mettre au point le programme et le cadrage. 


En fait nous gâchons de la mémoire car nous stockons dans un 
tableau une variable binaire : O ou 1. Nous n'utiliserions pas plus de 
mémoire en stockant dans chaque emplacement un nombre de 4 chif- 
fres dans lequel on pourrait lire 4 valeurs de la variable binaire. À vous 
de compléter le programme dans ce sens afin de profiter pleinement de 
la haute résolution de l'écran. 


Méthode 


Tout d'abord, nous considèrerons des équations r (z; a) telles qu'à 
une valeur de z corresponde une seule valeur der. Le cadrage est prévu 
pour un rayon maximal de 99 et une variation de z de O à 200. 


Recopiez donc le programme suivant en supprimant toutes les lignes 
dont le numéro ne se termine pas par zéro. 
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L'équation r (z, a) définit un volume que nous allons découper en 
tranches, donc à a constant. Nous décrivons les tranches en partant du 
bas (z = O) jusqu'en haut (z = 200). Nous commençons par les tranches 
de devant à gauche, pour aller derrière à gauche puis nous revenons 
devant à droite pour aller derrière à droite. Le solide est vu en perspec- 
tive, la rotation ah du regard revenant à faire un changement de variable 
dans l'équation : la tranche a devient la tranche 44=a—ah (ligne 20). 


Pour une valeur de z et de a on en déduit r donc les coordonnées 
«réelles» x, y, z du point (ligne 250). La rotation av du plan d'observa- 
tion modifie ces valeurs en y2 = y et z2 comme dans le programme 
«z (x, y}» (ligne 270). e est l'échelle de reproduction. z3 est l'ordonnée 
du point sur l'écran, 50 correspond au cadrage et peut être modifié 
suivant les besoins. Pour un essai, on n'utilise qu’un point sur 16 de 
l'écran, les coordonnées sont donc comptées de 4 en 4 (lignes 280 et 
310), y4 et z4 sont les indices du point sur l'écran. 


À une valeur de z correspond un point du bord de la tranche et un 
couple y4,z4 sur l'écran. Ces coordonnées sont mémorisées rp% et 
zp% (ligne 320). Si nous sommes sur la première tranche, le point 
est marqué et joint au précédent (ligne 340). En faisant varier z de 0 à 
200 nous marquons ainsi tout le bord de la tranche. Ce contour et 


2e tranche 


ZP% 


Fig. 1. 
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l'axe central définissent une surface fermée. En ligne 390 à 450 nous 
affirmons que cette tranche est opaque et que l’on ne pourra plus mar- 
quer de points à l’intérieur de cette surface fermée (figure 1). Dans 
cette première étude, rm% est constamment nul. Pour visualiser ce 
travail, vous pouvez ajouter la ligne 445 : 


PLOT 320 + 4.w.s, 4.t (figure 2) 


Break in 458 
Ready 
call 4106 


Fig. 2. 
L'ordinateur est en train de 
condamner les points de Ja 
seconde tranche. 





L'état de chaque point de l'écran (en fait un point sur 16) est noté 
p%ly4,z4). 

Lors du dessin de la tranche suivante on testera la valeur de p% 
avant de marquer un point (ligne 330). 


Les points à l'intérieur de la tranche sont répartis sur des rayons 
perpendiculaires dans la réalité à l'axe du solide, ils sont donc tracés 
à z constant ce qui nécessite un calcul assez pénible en lignes 420 à 450. 
Nous aurions pu balayer la tranche par des lignes horizontales sur 
l'écran (donc à z4 constant) mais la définition du début et de la fin de 
ces lignes est complexe. 


Le travail est fait pour la partie gauche de la figure (s = — 1) puis 
pour la partie droite (s = 1). 


Iustrations 


Après avoir vérifié le bon fonctionnement du programme et le 
cadrage, vous doublez la résolution avec : 


DIM p% (50,200) en ligne 130 
z4 = INT (zz/2) en ligne 310 


et en remplaçant 4z4 et 4 .(zz-z4) par 2.z4 et 2. (zz-z4) en ligne 340. 
Vous augmentez également le nombre de tranches en remplaçant 0,4 
par 0,2 ou 0,1 en ligne 200. 
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Nous vous proposons quelques équations en fin de programme, il 
suffit de les placer en lignes 40 à 49 pour observer la figure correspon- 
dante (nettoyez l’espace 40, 49 entre chaque essai). La «tour » nécessite 
un recadrage ; remplacez 50 par 150 en lignes 160 et 275. 


Figure 3 citrouille e=1,5 ah=0  av=45 


Break in 265 
Te (2 rs 

Ca a + = 

CZ, 

de TS 
Fig. 3. 
La moitié gauche de !a citrouille. 
(s=—1) 
tour e=0,8 ah=35 av=30 


= 


SR. 


LS 
Fig. 4. 
| La moitié droite de la citrouille. 
{s= 1) 
< 


L/ 


7h 
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Organisation du programme 


Lignes) 
40 équation de la surface 
60 introduction des données 
130 préparation des tableaux 


190 chaque côté 
200 chaque tranche 


250 où est le point ? 

330 est-il vu ? 

390 chaque tranche est opaque 
490 mise à zéro des p % 

540 exemples 


Liste des variables 


r,z,a,x,y coordonnées du point «dans la réalité » 

e échelle 

ah,av angles de vue 

sh,sv, cv lignes trigonométriques de ces angles 

aa angle a corrigé de ah 

sip%(y4,z4) = 1, il ne faut plus marquer de point à cet endroit 
u (z) ordonnée sur l'écran du point de l'axe central d’ordonnée z 
rp%lz) valeur maxmale de l’abscisse sur l'écran pour z et aa 
z2z3,k intermédiaires de calcul 

y4,z4 coordonnées du point sur l'écran 

yy,zz coordonnées du point précédent 

si d = 1 le point est vu 

w,t coordonnées du point écran rendu «opaque » 


Premier programme 


10 KEM'"surface r(z,a) “ 

20 REM 

30 CLS:GOTO 60 

40 r=(1-0.15#SIN(9#a) ) #SQGR(10000-(z2-100) 72) 


50 RETURN 

60 DEG 

70 PRINT: INPUT'"'e='";e 

80 PRINT: INPUT"ah (-90< <90 )='":ah 
90 PRINT: INPUT"av ( -90< <90 )=";av 
100 ah=ah+#*P1/180 

110 sv=SIN(av) :cv=C0S (av) 
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120 
150 
140 
145 
150 
160 
170 
180 
190 
£00 
210 
220 
225 
230 
240 
230 
260 
265 
270 
275 
280 
285 
290 
300 
310 
320 
330 


340 


CES: RAD 

DIM p4(50,100):DIM u(200) 
DIM rp4(200):DIM zp%4(200) 
DIM rmZ(200) 

FOR z=0 TO 200 
U(Zz)=INT((SO0O+e#z#cv)/4) 

NEXT z 

REM 

FOR s=-1 TO 1 STEP 2 

FOR aa=0 TO s*xPI STEP sx0.4 
a=aat+ah 

d=0 

c=0 

REM 

FOR z=Q TO 200 

GOSUR 40 
x=r#CO0S (aa) :y=rxSIN(aa) 
xb=r2+#C0S(aa):yb=r2+SIN (aa) 
Z2=2ZXCV-X#SV 
zb=50+ex(z+cv-xbxsv) 

Y4=INT (e*ARS (y) /4) 

by=INT (e*ABS (yb) /4) 
z3=50+e+*z2 

IF 2z3<0 OR z3>400 THEN GOTO 370 
z4=INT(z23/4) 
rpZ(z)=y4:zpl(z)=z4 

IF pZ(y4,z4)=1 THEN GOTO 560 
IF d=i THEN PLOT 520+4*%y4%s,4%x24: DRAUR 4%(yy-y 


4)#s,4#(z2z-24) 


350 
360 
362 
365 
364 


565 


d=i 

YyY=y4:2z=z4 

IF zb<O OR zb>400 THEN GOTO 370 
bz=INT(zb/4):rmX(z)=by 

IF pZ(by,bz)=1 THEN GOTO 367 

IF c=i THEN FLOT 320+4*%by+%s,4%bz:DRAWR 4%(ÿyp-b 


ÿY)#s,4#(zp-bz) 


366 
367 
370 
380 
390 
400 
410 
420 
430 
440 
450 
460 
470 
480 
490 
500 
510 
520 
530 
540 
550 
560 
570 
580 


c=i 

yp=by:zp=bz 

NEXT z 

REM 

FOR z=0 TO 200 

IF rpZ(z)=0 THEN GOTO 460 
k=(zp2(z)-u(z))/rp4(z) 
FOR w=rm£(z) TO rp2(z) 
t=u(z)+w#k 

pz(w,t)=1 

NEXT w 

NEXT z 

REM 

NEXT aa 

FOR w=0 TO 590:FOR t=0 TO 100 
pA(w,t)=0 

NEXT t:NEXT w 

[NEXT s 

END 

REM "tour 
b=100+70#CD0S(2#a) 

IF zùèb THEN r=80 ELSE r=200 
REM “champignon 

IF 2€80 THEN r=30:GOTO S0 
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590 IF z2100 THEN r=SGR(10000-—(z-100) 72): GOTO 50 


600 r=3.5#z-250 
610 REM "citrouille" 
629 r=(1-0.15#SIN(9#a) ) #SQR(10000-(z-100)"2) 


620 REM “tore" 

635 REM “tore" 

645 rr=SGR(2500-(z-100) 72) 
655 r=140+rr:r2=140-rr 


Améliorons le programme. 


De même que vous avez doublé le programme «z (x, y)» pour 
observer le dessous, de même en complétant votre programme r (z, a) 
par les lignes manquantes vous allez pouvoir observer l'intérieur de la 
surface r (z, a). Vous pourrez envisager alors les équations de surface 
telles que, pour une valeur de z et de a, il y ait deux valeurs possibles 


is 


Fig. 5. 










En figure 5, l'exemple du tore vous montre l'effet obtenu. La 
figure 6 précise les notations et voici enfin la correspondance, pas 
toujours heureuse il est vrai, entre les variables relatives aux deux parties 
de la surface. 


(en ligne 240 vous avez intérêt à faire varier z de 50 à 150 seulement.) 


EXTERIEUR  rp% 
INTERIEUR  rm% 


y4 
by 


z3 
zb 


z4 
bz 


r , 
r2 


x 
xb 


zZz 
zP 


ÿ 
yb 


yy 
yp 
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100 E — — — EE 


8 


50 
PR  rm% 


Fig. 6. 


Empilement de disques 


En permutant les variations de z et de aa, nous découpons le volume 
en tranches perpendiculaires à l'axe des z. Nous avons préféré vous 
donner le LIST complet bien qu'il n’y ait que très peu de modifications 
par rapport au programme précédent. Nous reprenons en figure 7 et 8 
les exemples déjà traités. 





Fig. 7. 
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10 
20 
50 


40 


50 
60 
70 


99 

100 
105 
110 
120 
1350 
140 
150 
160 
170 
180 
190 
195 
200 
219 
230 
250 
260 
270 


160 





Fig. 8. 


Deuxième programme 


REM‘"sur face r(z,a) en disques" 

REM 

CLS: GOTO 60 

r=(1-0.15#SIN(9#a)) #SQR(10000-—(2-100) 2) 


RETURN 
DEG 
PRINT: INPUT'"e=";e 
PRINT:INPUT'"ah (-90< <90 )=";ah 
PRINT:INPUT"av (0 < <90 )=";av 
ah=ah*P1/180 

g=SGN(av):1F g=0 THEN g=i 
sv=SIN(av):cv=C0S (av) 

CLS: RAD 

DIM pZ(50,100):DIM u(2Z00) 

DIM rp4(32):DIM zp4(32) 

FOR z=0 TO 200 STEP 5 
u(z)=INT((S0+e#z#cv)/4) 

NEXT z 

REM 

FOR s=-1 TO 1 STEF 2 

FOR z=100+g#100 TO 100-g#100 STEP -5%g 
d=0: FOR n=0 TO 31 
a=(n/19)+ah 

REM 

GOSUB 40 
x=r*COS(a):y=raASIN(a) 
Z2=ZA#CV-X#SV 


280 
290 
300 
310 
320 
330 


340 


YA4=INT (e*AERS(y)/4) 

z3=50+e*#z2 

IF 2z3<0 OR 23>400 THEN GOTO 3570 

z4=INT(z3/4) 

rpl(tn)=y4:zp4(n)=24 

IF pX(y4,z4)=1 THEN GOTO 360 

IF d=1 THEN FLOT 320+4+%y4%s,4#z4: DRAWR 4%(yy-y 


4)*xs,4#x(zz-24) 


250 
360 
370 
280 
390 
400 
410 
420 

30 
440 
450 
460 
470 
480 
490 
900 
510 
520 
330 
540 
5930 
560 
570 
580 
3970 


600 
619 
620 


d=i 

YY=Y4:zz=z4 

NEXT n 

REM 

FOR n=0 TO 31 

IF rpZ(n)=0 THEN GOTO 460 
k=(zpZ(n)-u(z))/rpz{(n) 

FOR w=0Q TO rpZ(n) 
t=u(z)+w#k 

pr(w,t)=i 

NEXT w 

NEXT n 

REM 

NEXT z 

FOR w=0 TO SO0:FOR t=0 TO 100 
pr(w,t)=0 

NEXT t:NEXT w 

NEXT s 

END 

REM "tour 

b=100+70#C0S (2#a) 

IF z>b THEN r=80 ELSE r=200 
REM “champignon 

IF z<80 THEN r=2z0:G0OTO 50 
IF z2>3100 THEN r=SORK(10000-(z-100)2):GOTO 50 


r=3.5#z-250 
REM “citrouille 
r= (1-0. 15#SIN(9%#a) )#SGKR(10000—(z-100) 2) 
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Représentation «Fil de fer» 
de polyèdres 


‘ Nous allons apprendre sur un exemple à représenter les arêtes 
de polyèdres convexes vus en perspective en supprimant les arêtes 
cachées. 


L'exemple choisi est un polyèdre à 10 sommets et 9 faces. Chaque 
face est un polygone plan. Ce polygone a au plus 4 sommets. Après 
avoir étudié cet exemple, vous pourrez réaliser sans peine des dessins 
de polyèdres à plus grand nombre de sommets, de faces, et dont les 
polygones ont un plus grand nombre de sommets. 


Définition du polyèdre 


En ligne 40 on prépare les tableaux de données : 10 correspondant 
au nombre de sommets, 3 au nombre de coordonnées (x, y, z), 9 au 
nombre de faces, 4 au nombre maximal de sommets pour une face. 


La figure 1 vous montre le polyèdre choisi qui représente une maison 
dans le repère Oxyz. Les numéros des sommets sont encadrés. 


En ligne 50 sont stockées les coordonnées des sommets dans l'ordre 
des numéros des sommets. Les groupes de 3 nombres sont séparés par 
une barre afin de rendre la lecture plus facile. Ces coordonnées sont 
introduites dans un tableau en ligne 80 . p (k, n) est la nième coordonnée 
du sommet numéro K. 


Les lignes 100 et 110 recadrent le polyèdre. Les valeurs 100 et 60 
peuvent être modifiées suivant les besoins. 
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40,140,120 





80,200,0 


Fig. 1. 


En ligne 140 sont stockées les définitions des faces. Par exemple, 
la troisième face est limitée par les sommets 2-3 7 6, l’ordre de ces nom- 
bres correspond à une rotation autour de la face. Si une vis était vissée 
dans cette face et si on la tournait dans ce sens, elle sortirait de la face 
vers l'extérieur du polyèdre (figure 2). Il en est de même pour les autres 
faces. Dans le cas des faces triangulaires deux sommets, les deux derniers 
du groupe de 4, sont confondus. En ligne 170, ces données sont intro- 
duites dans un tableau : f (kk,n) est le n ième sommet de la face Kk. 


(] 


el 
E G) 


Fig. 2. 
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Perspective 


En ligne 220 à 240 sont introduits les angles de vue, on calcule 
ensuite leurs lignes trigonométriques. En 270 sont calculés les coeffi- 
cients directeurs (non normalisés) de la direction R du regard (figure 3). 





Fig. 3. 


En ligne 240 est demandée l'échelle de reproduction que vous adap- 
terez à vos besoins. 


Regarder le polyèdre dans la direction R, c'est le placer dans un 
nouveau repère. Les coordonnées q (k,n) des sommets dans ce nouveau 
repère sont obtenues en ligne 300 à 320 par une matrice de rotation. Le 
travail est semblable à celui réalisé dans le programme «surface z (x, y)», 
nous n'y revenons donc pas. 


Arêtes cachées 


Soient À,, À,, À, À, les sommets d'une face. Nous allons mainte- 
nant calculer les coefficients directeurs de la normale à cette face. 
Cette normale est de même direction et de même sens que le produit 
vectoriel N = A, A, A A, A, (figure 4). Les composantes des vecteurs 
A; À, et A, A, sont calculées en lignes 370 à 420 et les composantes 
du produit vectoriel en lignes 430 à 450. 


Calculons enfin le produit scalaire T de ce vecteur N et du vecteur R 
(dans la direction du regard). En observant le signe de ce produit scalaire 
nous pouvons savoir si l'angle entre ces deux vecteurs est aigu ou obtus. 
S'il est aigu, nous voyons la face et nous la dessinons (lignes 490 à 510). 
S'il est obtus, c’est que la face est de l’autre côté du polyèdre, elle est 
cachée, nous ne la dessinons pas (ligne 480). Remarquons que dans ces 
calculs le module des vecteurs n'intervient pas, seuls interviennent la 
direction et le sens. 
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z} 


Aa A3 


Ai A2 


z} 
<{ 







angle aigu 
face vue 


& 





angle obtus 
face cachée 


Fig. 4. 


Z} 


illustrations 


Souvenez-vous enfin que l'imprimante déforme les figures (voir 
«graphes de fonction, … en polaires). 


Les illustrations obtenues correspondent à : 
Figure 5: ah=0 av=0 e=2 


6 : 90 0 2 
7 0 90 2 
8 45 0 2 
9 —45 —30 2 
10 30 30 1,8 
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Fig. 5. 


Fig. 6. 


Fig. 7. 
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Fig. 8. 


Fig. 9. 


Fig. 10. 
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Programme 


19 REM"polyedre" 
29 REM 

30 CLS 

40 DIM p(10,3):DIM q(10,3):DIM #(9,4) 

50 DATA 0,0,0,/,80,0,0,/,80,200,0,/,0,200,0,/,0,0, 
80,/,80,0,80,/,80,200,80,/,0,200,80,/,40,0,120,/,4 
0,140,120,/ 

60 REM 

79 FOR k=1i TO 10 

890 READ p{k,1),ptk,2),ptk,S),x# 

90 RENM 

100 p(k,2)=p{k,2)-100 

110 ptk,3)=p{k,3)-60 


120 NEXT k 

130 REM 

140 DATA 1,4,2,2,/,1,206:S53/523237563/535438373/34 
115558 /3 50675 9/75810,10,/,9,6,7,10,/,10,8,5,9 


, 
159 REM 
160 FOR kk=1 TO 9 
170 READ f(kk,1),f(kk,2),fCkKk,3),fCkKK,4),X# 
180 f(kKk,0)=f(kKk,4) 
190 NEXT kk 
200 REM 
210 DEG 
220 INPUT “ah=";ah 
30 INPUT "av=";av 
240 INPUT ‘"e='";e 
230 Ssh=SIN(ah):ch=C0S (ah) 
260 sv=SIN(av):cv=C0S (av) 
270 rx=cv#chiry=cv#shirz=sv 
280 REM 
290 FOR k=i TO 10 
300 x=p(k,1)#ch+p(k,2)#sh 
310 q(k,1)=-ptk,1)#sh+p{(k,2)#ch 
320 q(k,2)=-x#sv+p(k,3S)#cv 
330 NEXT k 
340 REM 
350 CLS 
360 FOR kk=1i TO 9 
370 x2=p(f(kk,2),1)-p{f(kk,1),1) 
380 x4=p(f(kk,4),1)-pCf(kk,1),1) 
390 y2=p{(f(kk,2),2)-ptf(kk,1),2) 
400 y4=p(f(kk,4),2)-pf(kk,1),2 
410 z2=p(f(kk,2),3)-p(f(kk,1),3) 
420 z4=p(f(kk,4),3)-p{f(kk,1),3) 
430 nx=y2#z4-22#y4 
440 ny=z2#x4-x2#z4 
450 nz=x2*+y4-y2#x4 
460 t=rx#nx+ry#Ny+rzHNz 
470 REM 
480 IF t<O THEN GOTO 530 
490 FOR n=0 TO 5 
900 MOVE 320+e#q{(f(kk,n),1),200+exq{(f(kk,n),2) 
510 DRAUR ex(q(f(kk,n+1),1)-q(f(kk,n),1)) ,e#(q(f(k 
kçn+1),2)-q{(f(kk,n),2)) 
S20 NEXT n 
530 NEXT kk 
540 END 
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Fractal 


Ce mot recouvre une technique permettant de dessiner des courbes 
et des surfaces presque émouvantes car elles imitent des reliefs naturels 
de façon étonnamment réaliste. 


Ces courbes ont un point commun avec les horizons ou les rivages, 
c'est que lorsqu'on les détaille (lorsqu'on augmente le nombre de points 
calculés sur une courbe allant d’un point donné à un autre point donné) 
la longueur de la courbe augmente et peut tendre vers l'infini. A la limite 
ces courbes peuvent remplir des surfaces (flocon de Von Koch, courbe 
de Péano). Le mathématicien Mandelbrot leur a attribué une dimension 
supérieure à 1 (une «fraction» d'où «fractal»). 


Remarquons tout d’abord que la fonction RND fournit à chaque 
mise sous tension de l’ Amstrad toujours la même suite de valeurs. Vous 
pouvez au besoin la faire tourner par : 


FOR X = TO. : Y = RND : NEXT X 


Courbe fractale 


Nous allons déterminer 257 points de cette courbe de x = 0 à 
x = 256. Le cadrage est prévu pour y variant de O à 400. || se peut 
toutefois que la courbe sorte de l'écran. Les valeurs 0, 200 en ligne 40 
et 99 en ligne 80 ne sont que des conseils. 


On demande d'abord de fixer les extrémités et le milieu de la courbe 
(ligne 110), ensuite on calcule le milieu de chacun de ces deux segments 
en ligne 150. Mais on s’écarte de ces points d’une quantité prise au 
hasard entre — e et + e (fig. 1).-On obtient ainsi deux nouveaux points 
de la courbe. On reprend le même processus sur les quatre segments 
obtenus, mais la quantité prise au hasard est comprise cette fois-ci entre 
— e/2 et +e/2, ainsi de suite. 
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abscisse  {n-1).x n.x {n+1).x 
valeur de t L 2 3 


Fig. 1 


Plus le segment à partager est petit, plus la quantité prise au hasard 
est petite. On peut parler d'échelle pour le hasärd. La courbe finale 
s'écarte d'autant plus de la courbe initiale que e est grand. 


S est le numéro de la boucle de calcul 
X est l'écart entre deux points 
T numérote les points 


Les exemples suivants ont été pris avec : 


Figure 2 : a=100 b=100 c=100 e=20 
3 100 100 100 90 
4 0 200 0 99 


Si la figure 2 évoque un horizon du massif central, la figure 5 évo- 
que plutôt les Alpes. 


Que TR NN RS “Éd CS 


Fig. 2. 


Fig. 3. 
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Fig. 4. 


Premier programme 


10 REM'"courbe fractale" 

20 REM 

359 CLS 

40 PRINT:PRINT'"nombres compris entre O0 et 200" 
90 INFUT'"'a="3;a 

60 INPUT"b=";b 

70 INPUT"c="3;c 

80 PRINT:INFUT'"degre (entre 0 et 99)";e 
90 CLS:PRINT'"calcul en cours" 

100 DIM y (256) 

110 y(0)=a:y(128)=b: y (256) =c 

120 t=%:x=64 

130 FOR s=1 TO 7 

140 FOR n=1 TO t STEP 2 

150 y(n#x)=(yC(n—-1)#x)+y((n+1)#x))/2 
160 yY(n#x)=y (n#x) +ex (—-1+2*RND) 

170 NEXT n 

180 t=2#*t+1:x=x/2:e=e/2 

190 NEXT s 

200 CLS 

219 FOR x=1i TO 256 

220 PLOT 80+2#x,y(x) +100 

230 DRAWR -2,y(x—-1)-y(x) 

240 NEXT x 

250 END 


Surface fractale 


Le «%» est destiné à alléger la mémoire de l'ordinateur. 


On représente la surface y (u, v),u et v étant compris entre O et 64. 
64 fois 64 valeurs de y sont donc à calculer. On demande d’abord à 
l'utilisateur de choisir les valeurs de y pour les sommets et le centre du 
carré u, v. Les autres valeurs sont calculées ensuite de la même façon 
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que pour la courbe. On détermine d’abord le milieu des droites côtés de 
chaque maille carrée puis leur centre comme le montre la figure appa- 
raissant sur l'écran. 

La surface obtenue est représentée en perspective suivant la techni- 
que évoquée dans le programme «représentation d’une surface z (x, y)». 
Nous n'y revenons donc pas ici. 


Nous avons fixé ah à 0. Nous n'avons pas envisagé de représenter le 
dessous de la surface mais vous savez comment le faire. 

Dans l’idée d'évoquer des paysages, nous avons choisi d'annuler tou- 
tes les valeurs de y négatives, ceci est fait en lignes 415 à 417, nous 
obtenons ainsi des images de rochers au bord de la mer. 


Les figures suivantes ont été obtenues avec : 


Figure 5 : a=—50 b=—-50 c=200 d=-50 e=-50 ee-80 av=20 
6 —50 —50 300 — 50 —50 200 20 
7 0 —100 0 200 —100 200 20 


Organisation du programme 


Ligne(s) 
10 introduction des données 
40 à 370 calcul des valeurs de y 
210 détermination des milieux des côtés 
310 détermination des centres des carrés 
390 tracé de la surface en perspective 


Remarque : dans ce deuxième programme, la valeur maximale de la 
quantité prise au hasard est notée ee. 
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Fig. 7. 


Deuxième programme 


10 REM'"surface fractale" 


20 REM 
30 CLS 

50 PRINT: PRINT"d e":PRINT:PRINT"  C':PRINT:FR 
INT'"a b" 

60 PRINT: INPUT"a="; a 

70 INPUT"b=";b 

80 INPUT"c=";c 

90 INPUT"d=";d 


100 INPFUT'"'e=";e 
110 PRINT: INPUT'"degre "see 
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120 DEG: PRINT: INPUT “angle (O0< <90)";av 

130 CLS 

140 DIM yZ(64,64) 

150 y2(0,0)=a:y2(64,0)=b 

160 yY2(0,64)=d:y2(64,64)=e 

170 t=1:x=52:PLOT 192,192 

180 FOR s=1 TO 6 

190 FOR v=0 TO 64 STEP 2+%x 

200 FOR n=1 TO t STEP 2 

210 y£(n#x,v)=((yZC(n-1)#x,v)+y2((n+1)%#x,vV))/2 )+e 
e#x(—-1+2#RND) 

220 yA(v,næ#x)=((yZCv, (n—-1)#x)+ÿYZ2(V, (n+1)#x))/7 )+e 
e#(—-1+2#RND) 

2350 PLOT 6#n#x,6#%v:PLOT 6#v,6#%n%x 

260 NEXT n 

270 NEXT v 

280 FOR n=1 TO t STEP 2 

290 FOR nn=i TO t STEF 2 

300 y4(n#x,nn#x)=(0.25#(y%4(nxx+x ,nnEx)+ÿ% (nHx-x , nn 


#x)+yY2(n#x ,NNnEX+X) +Y 2 (nEx ,NNnEx-x)))+ee*#(-1+2#RND) 
320 PLOT 6#n#x,6#nnxx 

330 NEXT nn 

540 NEXT n 


35350 y2(32,52)=c 
560 t=2#t+1:x=x/2:ee=ee/2 
570 NEXT s 
3580 REM 
590 CLS:FRINT"calcul en cours" 
419 Ssh=SIN(av):ch=CO0S (av) 
415 FOR w=0 TO 64:FOR z=0 TO 64 
416 IF y4(w,z)<0O THEN yZ(w,z)=0 
417 NEXT Zz:NEXT w 
20 CLS 
4:50 DIM m(65) 
440 FOR w=0 TO 64 
450 k=OQ 
460 FOR z=0 TO 64 
470 xe=6b*z+128 
480 ye=6#xsh#w+yZ(z,w)#ch 
490 IF ye<m(z) THEN GOTO 540 
9300 PLOT xe,ye 
910 m(z)=ÿe 
520 1F k=1 THEN DRAWR -6,yy-ye 
9330 k=1i 
940 yy=ye 
5930 NEXT z 
560 NEXT w 
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Méthode des moindres carrés 
Applications aux équations non linéaires 


En physique, les graphiques expérimentaux présentent en général 
des irrégularités. Si les points représentés se répartissent suivant une 
configuration assez régulière, le physicien ajustera la ou les grandeurs 
étudiées ; c’est-à-dire qu’il substituera aux grandeurs observées des 
grandeurs calculées. 

La première étape de l'ajustement porte sur le choix de la courbe 
d'ajustement, la seconde étape porte sur la détermination de l'équation 
de la courbe d'ajustement. Nous nous intéresserons aux équations de la 
forme : Y = f (X). 


Méthode des moindres carrés 


C'est une méthode analytique qui permet de calculer les différents 
paramètres W/ d'une équation de la forme Y = f (X). 

Représentons graphiquement les points expérimentaux de coordon- 
nées : 


(X1, Y1) IX2EN2E _ us (Xn, Yn) 
et les points de coordonnées qui se trouvent sur la courbe d'ajustement, 
soit (X1,Y'1) IX2V Z) ina (Xn, Y'n) 
Y v!, 
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La fonction Y = f (X) dépendant de plusieurs paramètres Wj 
des valeurs différentes de ces paramètres conduisent à des fonctions 
Y = f (X) différentes. 


La méthode des moindres carrés consiste à retenir la fonction 
Y = f(X) qui rend la plus faible possible la somme des carrés des écarts : 
Yi V'i. 

Il est hors de question pour nous de donner, dans cet ouvrage, la 


démonstration mathématique rigoureuse du calcul des différents para- 
mètres Wi de la fonction Y = f (X). 


Nous utiliserons la méthode de DEMING (1), reprise par WENT- 
WORTH (2), le lecteur intéressé pourra éventuellement consulter les 
articles des auteurs. 


Supposons que la fonction d'ajustement contienne la variable X et 
un certain nombre de paramètres : W1, W2 … WJ : soit Y’ (X, W1 … 
WJ). 


Il s'agit de trouver à l’aide d’un calcul itératif les valeurs optimisées 
des différents paramètres W/ (somme des carrés des écarts Yi — Y'; 
minimale). 


Soit AW/ la différence observée entre deux itérations successives 
pour le paramètre Wi. 


Ecrivons que : A W= Wi— Wi et Y'i= Y (Xi, Wi, W? .… Wj). 
à (Yi Yi) 


Posons F (W;) = 
ns F (Wi) | sw 


| SANT NE 


DEMING (1) a démontré que si les X/ sont connus avec précision, 
les A Wi se calculent en résolvant le système d'équations suivant. 


En posant Fi= Yi Y'i: 


J J J J 5 
Z F(WI.F(WT. Wi+ Z F(W1).F(W2) W2+..+ Z F(WI.F(WJ). WJ= 2 F(W1).F; 
1 1 1 


J J J o 
Z F(W2).F(W1). W1+ Z F(W2).F(W2). W2+..+ Z F(W2).F(WJ). WJ = È F(W2).F; 
1 1 1 1 


J J J J 
Z FIWJLF(W1). W1+ Z FIWJ).F(W2). W2+...+ Z F(WJ).F(WJ). W= ZX FIWJI.FS 
1 1 1 1 
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Liste des variables 


N nombre de points expérimentaux 

J nombre de coefficients de l'équation d'ajustement 
X (N) tableau pour les différents X; 

Y (N) LL LL LZA L22 Yi 


V(K) tableau pour les dérivés F (W) 

Remarque : V (1) = Y (1) — Y” 
Z (J, J) tableau pour les différents éléments de la matrice principale 
COR coefficient de corrélation 


Déroulement du programme 


Lignes) 


1 à 20 entrée des données (N et J) 

40 à 130 entrée de chaque couple de points expérimentaux (Xi, Yi) 
Remarque : à la ligne 70, on rentre les différentes valeurs 
des coefficients W; de départ (choix de l'utilisateur) 

140 et 150 on repère les valeurs maximales en valeur absolue des X; 
et des Y; 

200 à 250 initialisation des différents éléments de la matrice 

340 on pose V (1) = Y (1) — Y’ 
370 on pose DEF FN2Z (W) = Y (1) — Y’ 

390 à 410 les différentes dérivées F (Wi) = V (L + 1) sont calculées 
numériquement 


341 on calcule le coefficient de corrélation COR après cha- 
que itération 
500 saut au sous-programme d'inversion de matrice 


Les calculs se déroulent sous vos yeux jusqu'à la convergence du 
coefficient de corrélation (noter les différents W/ après avoir arrêté le 
programme en appuyant 2 fois sur ESC à la convergence). 


1000 à 1200 l'ordinateur efface l'écran. Il demande à nouveau les 
différents W/ calculés à la convergence 
1025 on trace les axes, les points expérimentaux à l'échelle 
et la courbe d'ajustement 
Remarque : en ligne 1080, n'oubliez pas de rentrer la 
fonction d'ajustement Y” 
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1 CLS 

S REM MOINDRES CARRES 

6 REM 4636262666 DE EE 

7 REM 

10 INPUT "NB DE POINTS EXP.N=";N 

20 INPUT "NB DE COEFFTS.J=";7J 

35 DIM X(N):DIM Y(N):DIM Z(J,J):DIM D(J) 
40 M=0:P=0: DIM W(J+1):DIM R(J) 

42 K=J+1:DIM V(K) 

44 DIM U(J,J):DIM E(J):DIM F(J) 

930 DW=0.001:PR=O,.00000001:MI=1E-50 

99 MODE 2 

60 FOR 1=i TO J 

70 INFUT "W="3;WCI): PRINT "W's;1l;"="; WI) 
80 NEXT I 

90 CLS 

100 FOR I1=1i TON 

110 PRINT "I="3s1 

120 INPUT "X=";X(C1): INPUT "Y="3; Y(I) 

130 PRINT "X";l;"="sX(IDg "ets" Ytals "=" Y (I): FRINT 
140 IF ABS(Y(I1))>P THEN F<ARS(Y(1)) 

150 IF ABS(X(1))2M THEN M=ABS(X(I)) 

160 NEXT I 

190 S=0 

200 FOR C=1 TO J 

210 FOR L=1i TO J 

220 Z(L,C)=0 

230 NEXT L 

240 NEXT C 

230 FOR L=1i TO J:D(L)=0:NEXT L 

300 RKEM CALCUL ITERATIF 

3510 REM 363626266666 6 NE 

316 FOR I=1 TO N 

3209 X=X(I) 

340 V(1)=Y(I1)-(W(1)#(1-X)-W(1)#(1-X)#EXPF(-W(ZD#X/( 
1-X))+W(3)#X-W(3S)#X#EXF(-W(4)#(1—-X) /X)) 
341 D=V(1)#V(1):S=S+D 

3350 FOR L=i TO J 

360 W=W(L) :A=W(L) 

370 DEF FNZ(W)=Y(I)-(W(1)#(1-X)-W(C1)#(1—-X)#EXF (—W( 
2) #X/(1-X))+W(3)#X-W(5)#X#EXF(-W(4)#(1—-X) /X)) 
390 W(L)=A+DW: Zi=FN Z (W) 

400 W(L)=A-DW: ZZ=FN Z(W) 

410 V(L+1)=(Z1-22) / (2#DW) 

411 B=V(L+1)#V(1) 

412 D(L)=D(L)+B 

415 W(L)=A 

420 NEXT L 

440 FOR C=1 TO J 

445 FOR L=1i TO J 

450 H=V(L+1)#V(C+1) 

4355 Z(L,C)=Z(L,C)+H 

469 NEXT L 

470 NEXT C 

490 NEXT I 

491 CR=S/N 

492 FRINT “"COR=";SQR(CR) 

4953 PRINT 

300 GOSUB 2021 

600 FOR L=1 TO J:R(L)=0 
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610 F 
629 KR 
630 N 


OR C=1 TO J 
(L)=R(L)+Z(L,C)#D(C) 
EXT C 


649 NEXT L 


641 F 
642 
6453 N 


OR F=i TO J 
(F)=W(F)-R(F)2:PRINT W(F) 
EXT F 


660 GOTO 190 


1000 
1010 
1020 
1025 
10370 
1049 
1041 
1942 
1043 
1050 
1060 
1070 
1080 
3) #X- 
1085 
1090 
1100 
1200 
2000 
2010 
2021 
2022 
2023 
2925 
2026 
2030 
2031 
2035 
2956 
29037 
2040 
2045 
2046 
2050 
2060 
2070 
2080 
2100 
2110 
2120 
2150 
2140 
2150 
2160 
2170 
2180 
2199 
2200 
2210 
2220 
2230 
2249 


CLS 

FLOT 0,200:DRKAWR 620,0 

PLOT 0,0:DKAWR 0,400 

FOR F=1 TO J:INPUT W{(F):NEXT F 

FOR 1=1i TO N 

PLOT X(1)/M#600,Y(1)/F#180+200:DRAWR 2,0 
PLOT X(1)/M#600,Y(1)/P#180+200:DRAWR 0,2 
PLOT X(1)/M#600,Y(1)/P#180+200:DRAWKR -2,0 
FLOT X(1)/M#600,Y(1)/P#180+200:DRAWR O,-2 
NEXT I 

FOR K=i TO 600 

X=(M/600) #F. 
Y=W(1)#(1-X)-W(1)#(1-X)#EXF(-W(2)#X/ (1—-X) )+W( 
W(3)#X#EXF(-W(4)#(1—-X)/X) 

IF ABS(Y/P#180)2>200 THEN GOTO 1100 
FLOT K,Y/P#180+200 

NEXT K 

STOP 

REM SP, INVERSION DE MATRICE 

REM 63636626 DD DEEE ED DEEE DD DE 

FOR C=1i TO J 

FOR L=1i TO J 

U(L,C)=Z(L,C) 

NEXT L 

NEXT C 

DET=1 

FOR G=1 TO J 

MA=MI 

FOR L=G TO J 

FOR C=6G TO J 

A=Z (L,C) 

IF ABS(A)<=MA THEN GOTO 2050 

IM=L : JM=C: MA=ABS (A) 

NEXT C 

NEXT L 

IF MA<2>MI THEN GOTO 2100 

STOP 

PT=Z(IM,JM) : DET=DET#PT 
F(G)=1IM:E(G)=JM 

IF IM=G THEN GOTO 2160 

FOR C=1 TO J 
A=Z(IM,C):Z(IM,C)=2(G,C):2(G,C)=A 
NEXT C 

IF JM=G THEN GOTO 2200 

FOR L=i TO J 
A=Z(L,G):Z(L,G)=2Z(L,JM:Z(L,JM)=A 
NEXT L 

FOR L=1 TO J 

IF L=G THEN GOTO 2300 

R=Z (L,G)/FPT 

FOR C=1 TO J 

IF C=G THEN GOTO 2290 
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2250 DD=Z(L,C)-R#Z(G,C) 
2260 1F AES(DD)<ABS(Z(L,C))#FR THEN DD=0 
2270 Z(L,C)=DD 

2290 NEXT C 

2300 NEXT L 

2400 FOR Z=1 TO J 

2410 Z(Z,G)=Z(Z,G)/FT 

2420 Z(G,Z)=-2(G,2)/FT 

2430 NEXT Z 

2500 Z(6,G)=1/PT 

2510 NEXT G 

2520 FOR G=1 TO J 

25350 Z=J-G+1 

2540 FL=F (2) 

2550 IF PL=Z THEN GOTO 2600 

2560 FOR L=i TO J 

2570 A=Z(L,2Z):2(L,2Z)=Z(L,PLI:Z(L,PL)=A 
2580 NEXT L 

2600 FC=E(Z) 

2610 IF PC=Z THEN GOTO 3000 

Z620 FOR C=1 TO J 

2620 A=Z(Z,C):2Z(Z,C)=Z(FC,C):Z(FC,C)=A 
2640 NEXT C 

3000 NEXT G 

3010 RETURN 


Nous donnons 3 exemples d'utilisation de ce programme. 


Exemple 1 : la fonction d'ajustement est une droite. 


kr 
.848528138 


Sun 
S- : 
rt 
sé 


Out 
=. - 
11 


0.848528138 


Q 


0.848528138 


Sur 
=. 


[e] 
Out 
=. 
11 œ 


0.848528138 


6.848528138 


Fo Se 
De - 


CEE 


in 41 


2% co 
Ex 


O% 
ES 
; 
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Y'=W{1)+X + W (2) 
Rentrer au départ W (1) = 1 et W{(2) = 1 





Exemple 2 : la fonction d'ajustement est une fonction du cinquième 
degré rendant compte de l’enthalpie de «mélange pour le propanol- 
eau. 


Mélange 1 PROPANOL 
Mélange 2 EAU 
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LEYDET (3) a montré que les fonctions de mélange en thermody- 
namique pouvaient être représentées par des fonctions d'exponentielle, 
Essayez les mêmes valeurs itérées avec la fonction donnée en lignes 340- 
370 et 1080 dans le listing. 


Remarque : les fonctions de mélange s’annulent pour x = 0 et x = 1. 





Y' = Xx(1—X)(W (1) +W(2)eX + W(3) # XT2 + W (4) # XT3 + W (5) + xf 4 
+ W1(6) « XT5) 


On part de W (1) = W1{(2) = .... = W(6) = 1 
Essayer ensuite la fonction : (donnée dans le listing) 


Y'=W{t1)e(1— X)— W (1) » (1 — X) « Exp (—W (2) » X/(1 — X)) 
+ W(3) »# X —W (3) * X » Exp (— W (4) » (1 — X}/X) 


partir de : W (1) = — 100 W (2) = 30 
W (3) = 300 et W (4) = 3 
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Exemple 3 : lissage d'une courbe de températures relevées à Marseille au 
mois de juillet. 


80.9757324 ® 
26: 3172363 x (heures) y ( C) 

3.78041077 6 16,5 
8 .202210903 a 

3.63742È-03 8 , 
ch 8.243686623 10 23,8 

12 28,4 

98.9628921 ; 
AU 16 + 
RCE 18 23,1 
co 0.142433226 20 5: 

89.8292149 22 : 
-26.9036656 

3.77878355 
6202121507 

3.635825 -0 


3 
COR= 6.343493402 


Break in 2268 
Ready 
CALL 44008 





Y'= WU) +W(2)» X + W(3)+ X12+W(4)» X13 +W(5)+X7T14 
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Voici des programmes de mathématiques écrits en BASIC. 


Accompagnés de leurs commentaires, ils transforment l'ordinateur 
en un outil pédagogique de valeur et un instrument de calcul puissant 
pour résoudre bien des problèmes. 


En second plan, les auteurs ont voulu montrer que les mathéma- 
tiques ne sont pas mortes avec l'usage de l'ordinateur, bien au contraire : 
l'écriture d'un programme comme la critique des résultats nécessite 
des connaissances mathématiques solides et variées. 


Cet ouvrage original devrait intéresser et surtout rendre service 
au lecteur, pourvu qu'il ait le niveau mathématique d'une classe termi- 


nale scientifique. 
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